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Dedicated to A.J. Coleman 

on the occasion of his 80-th birthday 

Preface 

Green's functions are widely used in various branches of mathematics and theoret- 
ical physics. Despite its rather long history, the method of Green's functions continues 
to develop. This development concerns first of all the application of this technique 
to statistical mechanics. The latter is, probably, the most complicated subject which 
the method of Green's functions has been applied to. Generally speaking, practically 
all parts of theoretical physics can be treated as subsections of statistical mechanics 
(this statement is explained more precisely in Appendix). It is just the complexity 
of statistical mechanics that perpetually forces us to develop the method of Green's 
functions. 

There exists a volumnous literature on this method and its various applications. In 
no way my aim was to give a review of all that. The aim of this work is to describe 
the basic mathematical properties of Green's functions used in statistical mechanics, 
as well as the equations defining these functions and the techniques of solving these 
equations. No particular applications are touched in order not to be distracted from 
the main aim of analysing the underlying mathematical structure. Such a concentrated 
analysis of the fundamentals of the method may, to my mind, help to grasp its most 
important features and to better understand its strong and weak points. This could 
also give hints of how to improve and develop the approach for overcoming existing 
difficulties. One such a development that has not yet been published in a review form 
is presented here in Part 3. This is Correlated Iteration Theory that differs from all 
other known variants of perturbation theory for Green's functions by the combination 
of two factors: the systematic formulation of an algorithm for obtaining subsequent 
approximations and the consistent consideration of interparticle correlations at each 
step of the procedure. 

The formalism of the correlated iteration theory is based on general properties of 
Green's functions, which are expounded in Part 2. Without studying these general 
properties it is impossible to proceed further. And, in turn, to introduce Green's func- 
tions for statistical systems and to study the properties of these functions, it is necessary 
to have those main definitions that make the foundation of statistical mechanics. These 
basic notions are formulated in Part 1. 

Statistical Green's functions are a natural generalization of correlation functions 
and reduced density matrices. The relation between these objects is repeatedly traced 
throughout the text. 
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1 Basic Notions 



This part contains those main definitions which statistical mechanics is based on and 
which are necessary for understanding the following parts. 

1.1 Pure States 

Consider a system of N objects which we shall traditionally call particles. Put into 
correspondence to each particle a variable x G V whose particular nature is, for the 
time being, of no importance. This can be a multivariable incorporating the real-space 
coordinates, spin, isospin and other degrees of freedom, such as colour and flavour. 

For dealing with N particles, we need N such multivariables. Enumerating the 
variables (but not the particles) by the index i = 1, 2, . . . , N, we have the set of iV 
variables 6 V,. For brevity, we shall use the notation 

x n :— {xi\ i — 1,2, . . . ,n < N} (1.1) 

for the set defined on the manifold 

F:=? 1 xV 2 x...xV„. (1.2) 

In addition to variables (1.1) we need one more variable 

teR + := [0,oo) 

called time. 

Let functions *ff(x N ,t), given on the domain x M + , pertain to a Hilbert space 
Hxf, with the scalar product 

(tf, #') : = J ^(x", t)*'(x N , t)dx N , (1.3) 

in which ^* means complex conjugate of integration goes over Y N and the differ- 
ential measure on (1.2) is denoted by 

dx n := dx\dx2 ■ ■ ■ dx n . 

As usual, the integration implies an actual integration over continuous manifolds and 
summation over discrete ones. 

If the elements of the Hilbert space Ti^ are normalized, 

= (1-4) 

then the functions \& are called wave functions or quantum states; respectively, the 
space Hq, is the space of quantum states. The quantity 

\*(x N ,t)\ 2 := V*(x N ,t)V(x N ,t) (1.5) 
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is interpreted as the distribution of particle density. 

The particles are called identical or indistinguishable if their corresponding variables 
pertain to the same manifold, Vj = V, and there is no difference which variable is 
ascribed to which particle, so that the particle-density distribution (1.5) does not 
change under the permutation of any two variables: 

iP(* 4 ,z^(^)rH^V)r, (i-6) 

where P(xi,Xj) is the operator permuting the variables x iy Xj G V. Emphasize that it 
only makes sense of permuting variables but not particles. 

The simplest sufficient condition for (1.6) is the symmetry principle 

Pfa, x j )^(x N , t) = ±^(x N , t), (1.7) 

stating that the quantum states associated with identical particles are either symmetric 
or antisymmetric with respect to the permutation of variables. Symmetric quantum 
states describe Bose particles or bosons while antisymmetric quantum states, Fermi 
particles or fermions. 

A set {A} of self-adjoint linear operators A(x N ) is given on the space of quantum 
states TCq,, such that the mathematical expectations 

(A):=(*,Atf) (1.8) 

define observable quantities. A linear algebra of these operators forms the algebra of 
observables {A}. The set {(^4)} of mathematical expectations (1.8) is called observable 
state. 

One of the operators from the algebra of observables plays a special role. This 
is a self-adjoint operator H whose mathematical expectation gives the energy of the 
system. The operator H, called Hamiltonian, defines the evolution equation 

ih^(x N , t) = H(x N )V(x N , t) (1.9) 

for quantum states. This is the Schrddinger equation, h being the Planck constant. 

Using wave functions one may construct another type of objects, called density 
matrices [1-3]. Introduce the pure density matrix 

(1.10) 
whose elements 

D N (x N ,y N ,t) := V(x N : t)V*(y N : t) (1.11) 
are products of wave functions satisfying the symmetry principle (1.7). Consequently, 

P(x l ,x J )D N (x N ,y N ,t) = ±D N (x N ,y N ,t). (1.12) 



D N (t) 



D 



N 



(x N ,y N ,t) 
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The trace of (1.10), in compliance with (1.4), is 

7VA,(t):=(tt,tt) = l. (1.13) 
Also, introduce the operator matrices 

A N = [A N (x N ,y N )] (1.14) 

with elements 

An(x n , y N ) ■— A(x N )S(x N — y N ), (1.15) 

where 

5{x n - y n ) := 5(x 1 - y 1 )6(x 2 - y 2 ) ■ ■ ■ S(x n - y n ) 

is a product of the Dirac 5-functions for continuos variables and of the Kroneker 5- 
symbols for discrete. Using (1.10) and (1.14), the mathematical expectation (1.8) is 
written as 

(A) — Tr A N D N (t). (1.16) 
The evolution equation for (1.10) follows from the Schrodinger equation (1.9), 

thj t D N (t) = [H N ,D N (t)}, (1.17) 

where [ , ] means the commutator 

[A, B] := AB-BA. 

Eq. (1.17) is the Liouville equation. According to (1.17), the evolution equation for 
(1.16) is 

r\ r)A 

Mj t (A)=ih(^) + ([A,H}). (1.18) 

The pure density matrix (1.10) and operators of observables (1.14) permit several 
representations. To pass from the ^-representation given by (1.11) and (1.15) to an- 
other representation, one chooses a complete basis {^/ q (x N )} of orthonormal quantum 
states in Hy, 

where q is a multiindex labelling quantum states. With an expansion 

*(x N ,t)=Y,C q (t)%(x N ), 
i 

Eq. (1.4) gives 

Eic*(*)l 2 = i- 

Instead of (1.10), one may present the pure density matrix as 

D N (t) = [D qp (t)\ (1.19) 
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with elements 

D qp (t) := C q (t)C* p (t), (1.20) 

which again have the form of a product, the function C q (t) playing the role of the wave 
function in the g-representation. In the mathematical expectation (1.16) one has to 
put now, instead of (1.14), 

A N = [A gp ] (1.21) 

with elements 

A qp :=(^ q ,A^ p ). (1.22) 

In any representation, the pure density matrix consists of elements that are factoriz- 
able as in (1.11) or (1.20). Because of this and of normalization (1.13), the convolution 
property 

J D N (x N , z N , t)D N (z N , y N , t)dz N = D N (x N , y N , t) (1.23) 
holds, which can be written in an operator form as 

D 2 N (t) = D N (t). (1.24) 

The latter defines the pure density matrix as an indempotent operator. 

A pure state is an observable state {(A)} consisting of mathematical expectations 
(1.16) defined through a pure density matrix. The indempotency property (1.24) is the 
necessary and sufficient condition of a pure state. 

Reducing the number of variables in (1.11) by integrating 

D n (x n , y n , t) := J D N (x n , z N ~ n , y n , z N ~ n , t)dz N ~ n , (1.25) 

one defines the reduced density matrices, or submatrices, 

D n (t) := [D n (x n , y n , t)} (n — 1,2, ... , N), (1.26) 

implying that Dq := 1. It is not important which of the variables are integrated out in 
(1.25) since (1.11) is symmetric under the double permutation 

P(x i ,x j )P(y i ,y j )D N (x N ,y N ,t) = D N (x N ,y N ,t). 

The normalization for (1.26), 

Tr D n (t) := j D n (x n , x n , t)dx n = 1, (1.27) 

follows from (1.13). The elements (1.25) do not satisfy the convolution property (1.23), 
except for the trivial cases n = or n = N. Thus, the reduced density matrices (1.26) 
do not represent indempotent operators, 

D 2 n (t)^D n (t) (n = l,2,...,N-l). 
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The reduced density matrices are useful for dealing with operators depending only 
on a part of the variables x N . For example, for an operator A(x N ) defining the operator 
matrix 

A n = [A n (x n ,y n )} (1.28) 

with elements 

A n (x n ,y n ) :=A(x n )5(x n -y n ), (1.29) 
one can write the expectation (1.16) as 

(A) = Tr A n D n {t). (1.30) 

In this way, to calculate (1.30) for an operator (1.28) whose elements (1.29) involve 
n < N variables, one needs to know only the corresponding reduced density matrix 
(1.26). Reduced density matrices have become a very convenient tool for the many- 
body problem [4-14]. 

1.2 Mixed States 

In reality, any considered system of iV particles is always a subsystem of a bigger 
system. Real observable systems can be never completely isolated [15-19]. Therefore, 
in addition to the variables x N ascribed to the subsystem of interest, there exists a set 
X of other variables associated with the rest of a bigger system, which can be called 
the background. The whole system consisting of a subsystem and a background can be 
named the universe. 

Assume that the observable state of the universe is a pure state characterized by a 
pure density matrix 

D(t) = [D(x N ,X,y N ,Y,t)}, (1.31) 

whose elements 

D(x N , X, y N , Y, t) := *(x N , X, t)^(y N , Y, t) (1.32) 

are products of the universe wave functions. 

The algebra of observables {A} related to the subsystem of interest consists of oper- 
ators A(x N ) depending on the subsystem variables x N . The mathematical expectations 
can be defined as in (1.16) but with the reduced density matrix D N (t), whose elements 
are 

D N (x N , y N , t) := J *(x N , X, t)**(y N , X, t)dX . (1.33) 

To pass to another representation, one has to choose an orthonormal basis {^ q (x N )} 
and to make an expansion 

y(x N ,X,t)=Y,C q (X,t)^ q (x N ). 

Then, (1.33) may be written as 

D N (x N ,y N ,t)=J2D qp m q (x N )%(y N ), (1.34) 

qp 
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where 



D qp (t) = J c q (x,t)c;(x,t)dx 



are the elements of the reduced density matrix D^{t) in the g-representation. 

There is a sense in considering a system of N particles as a subsystem of the universe 
solely if any of these particles has ever interacted with its surrounding background. This 
is just what always happens in reality. It would be an extremely unrealistic fantasy 
to imagine that a considered subsystem has never interacted with anything else. The 
fact that a subsystem cannot be completely separated from its background means that 
the universe wave function ty(x N , X, t) cannot be factorized into a product of functions 
each depending separately on x N and X. This is called the principle of nonseparability 



If the universe wave function is not factorizable in the above sense, then the elements 
of the reduced density matrix D N (t), given by (1.33), are not factorizable into a product 
of functions depending separately on x N and y N , as in (1.10). Therefore (1.33) cannot 
define a pure density matrix. 

A density matrix consisting of nonfactorizable elements will be called the mixed 
density matrix. The observable state {(A)}, formed by mathematical expectations (A) 
which are defined through a mixed density matrix, is named the mixed state. 

The necessity of characterizing any system of particles by a mixed state is not 
based only on the principle of nonseparability, as is described above, but can even be 
strengthen. Actually, the notion of an observable state itself presupposes the possibil- 
ity of measuring observable quantities. But the process of measurement is impossible 
without interactions between the considered system and a measuring apparatus. Con- 
sequently, the latter is to be treated as a part of the background from which the system 
cannot be isolated. 

Moreover, the concept of isolatedness of a subsystem from the rest of the universe is 
contradictory in its own, since to realize the isolation and to check that the subsystem 
is kept isolated, one again needs to use some technical and measuring devices [23,24]. 

Suppose that to check the isolatedness of a subsystem we use such a gentle procedure 
which involves only the minimally disturbing measurement [23,24], defined as follows: 
A measurement, finding the subsystem in a quantum state *& q (x N ), leads solely to the 
appearance of a phase factor e^ q , 



where £ 9 is a random real number, £ g e [0, 2n]. The measuring process (1.35) does not 
change the distribution of particle density ^(a^)) 2 in the quantum state q, as well as 
the matrix elements (^f q , A^f q ). However, the total wave function 



acquires the dependence on a set ( := of random variables, which plays here the 
role of the set X of background variables. To preserve the normalization 



[13,20-22]. 





*(x N ,Z,t)=J2C q (t)e^ q (x N ) 



(1.36) 




(1.37) 
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the differential measure dt; is to be defined as 

2tt ' 



^ : =Il!r> e 9 e[0,27r]. (1.38) 
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The elements of the reduced density matrix, according to (1.33), are given by 

D N (x N , y N , t) := J *(x N , £, t)^(y N , £, t)d£ . (1.39) 

The mathematical expectation (A) of an operator A(x N ) is Tr A^D^{t), where the 
elements of An and D^{t) are defined either by (1.15) and (1.39) or by (1.22) and 

D qp (t) = \C g (t)\ 2 6 qp , (1.40) 

respectively. The latter follows from (1.36), (1.38) and from the integral 



Neither (1.39) nor (1.40) are factorizable. Hence, the subsystem is in a mixed state 
even for such an ideal case of the minimally disturbing measurements. 

The above discussion forces us to accept the principle of nonisolatedness which tells 
that no system can be completely isolated and has to be described by a mixed density 
matrix. This general statement, however, does not forbid the possibility that for some 
particular cases a pure density matrix could be a reasonable approximation. 

We shall call a system quasi-isolated if it is characterized by a mixed density matrix 
and its time evolution is governed by a Hamiltonian H(x N ) depending on the system 
variables x N , but not on the background variables X, 

ih^(x N , X, t) = H(x N )V(x N , X, t). (1.41) 

The mixed density matrix of a quasi-isolated system satisfies the same Liouville equa- 
tion (1.17) as a pure density matrix. 

The formal solution of the Liouville equation may be written as 

D N (t) = U N (t)D N (0)U+(t), (1-42) 

where U^it) is the evolution operator which is a unitary operator obeying the 
Schrodinger equation 

ih^U N (t) = H N U N (t) (1.43) 
at 

and having the properties 

U N (t)U N (t') = U N (t + t'), 
U^(t) = U N 1 (t), U N (0) = 1. (1.44) 
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In a particular case, when does not depend explicitly on time, the evolution oper- 
ator, according to (1.43), is 

U N (t) = exp (~H N (t)) = o) . (1.45) 

With the help of the evolution operator the time evolution for the operators of 
observables is given as 

A N (t) := U+(t)A N U N (t). (1.46) 

The form A^(t) is named the Heisenberg representation for the operators, as compared 
to that of A N which is called the Schrddinger representation. From (1.46) and (1.43) 
the evolution equation 

d d 

ih-A N (t) = ih—A N (t) + [A N (t), H N (t)} (1.47) 

follows, in which 

g-A N (t):=U+(t)-^U N (t). (1.48) 

Eq.(1.47) is called the Heisenberg equation. Its important particular case is when the 
operator A N does not depend explicitly on time, then (1.47) and (1.48) yield 

ihj t A N (t) = [A N (t), H N (t)} = o) . (1.49) 

Using the properties of the evolution operator and the equality Tr AB = Tr BA, 
one can write the expectation value (1.16) in one of two forms 

(A) = Tr D N (t)A N = Tr D N (0)A N (t). (1.50) 

Respectively, the evolution equation for (1.50) may be presented in two ways 

ihj t (A) = i h(^)+Tr [H N ,D N (t)]A N = 
r) A 

= ih (^r) + Tr D N (0)[A N (t), H N (t)]. 
ot 

The identity of these two forms results from the properties of trace and from the relation 
(1.46) between the Schrodinger and Heisenberg representations. 

For mixed states, the reduced density matrices can be introduced as for pure states 
in (1.25) and (1.26), though now there is no principal difference between D N {t) and 
D n (i) with n < N, since for a mixed state D N (t) itself is a reduced density matrix. 
Dealing with mixed states it is reasonable to consider all density matrices D n (t) as 
reduced, with the index n ranging from zero, when D (t) := 1, up to infinity. 

The reduced density matrices for mixed states have the symmetry property 

P(^,^)L>„(x",y n ,0 = ±D n (x n ,y n ,t) 
as (1.12) for pure states. 
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1.3 Field Operators 



A powerful method of treating quasi-isolated systems is the field approach based on 
the notion of field operators. This approach is also known as the method of second 
quantization. There are several ways of describing the latter. Here we expound a 
variant which is, to our mind, the most general and best mathematically grounded 



Let 7i n be a Hilbert space of square-integrable functions h n (x n ) with a domain V" 
and the scalar product 

(K,g n ) := Jh* n (x n )g n (x n )dx n , 
where the integration goes over V n . Construct a Hilbert space H as a direct sum 



H := (B^Hn (Ho := C). 
A vector h G H is presented as a column 



h 
/ii(xi) 
h 2 (x 1 ,x 2 ) 



;i.5i) 



h := [h n {x n )\ = 



h n {X\) X2-, ■ ■ • , Xn) 



(1.52) 



The latter is a quantum state, or microstate, and 7i is a space of states or phase 
space. The scalar product in (1.51) is given by 

oo 

(h,g) := h + g = Y( h n,9n)- 

n=0 

Define matrices of the form 

A := [A mn {x m , y n )}, (1.53) 
whose action on the vector (1.52) is a vector 



Ah = 



oo „ 

]T / A mn (x m ,y n )h n (y n W 

,n=0 J 



and the multiplication of two such matrices is given by the rule 



AB 



J2 A mk (x m ,z k )B kn (z k ,y n )dz k 
.k=o J 
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Introduce the annihilation operator 



m n—1 



y/n5( 



x 



n—1 



y n-l )6{x _y n) 



ip(x) = 5. 
and the creation operator 

ft(x) = S mn+ Wn+ 15 (x n - y n )5(x n+1 - x) 



1.54) 



(1.55) 



having the form of the matrix (1.53) with the corresponding multiplication rules. As 
follows from above, the operators (1.54) and (1.55) are adjoint to each other; these are 
the field operators. Since the matrices (1.54) and (1.55) contain the Dirac ^-functions, 
the field operators are operator distributions. The correct definition of distribution 
presupposes the existence of a space of square-integrable functions f(x), so that 



i>f = / f(x)^(x)dx, 



f(x)ft(x)dx 



(1.56) 



becomes operators in 7i, whose domain of definition, as well as their range, contain the 
set of all finite vectors. 

Action of the field operators on the vector (1.52) gives 



tjj(x)h = Vn + lh n+ i(x n , x) 

$(x)h = \/nh n -\{x n ~ 1 )5{x n — x) 

so that tpfh and ip^h are really finite vectors. 

The explicit definition of the field operators (1.54) and (1.55) imposes no simple 
commutation relations on these matrices. To understand this, it is sufficient to examine 
their products 



ip{x)ip{x') 



S m n -2^Jn{n - l)5{x n - y n 2 )5{x - y n -i)5{x' - y n ) 



S m n+2 J{n + l)(n + 2)8{x n - y n )5{x n+l - x')5{x n+2 - x) 



^(x)^(x') = [S mn (n + l)5(x n - y n )5(x - x')} , 
^(x')^(x) = [5 mn n5{x n - 1 - y^Sixn - x')5(y n - x)} , 

The algebraic properties of operators depend on three factors: the explicit form of 
an operation, the properties of the space on which the operators are defined, and the 
existence of some additional conditions. This is well illustrated by the example of field 
operators. Their explicit definition as matrices (1.54) and (1.55) does not permit to 
extract any simple commutation relations. Even if we narrow the domain space to a 
subspace enjoying some simple symmetry properties, this does not necessary lead to 
simple commutation rules for operators. The latter can be demonstrated by the actions 



ip(x)ip(x')h 



^/{n+l){n + 2)h n+2 {y n , x, x') 
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ft(x)ft(x')h = ^Jn(n-l)h n ^ 2 (x n 2 )5(x n ^ - x')5(x n - x) 

r/>(x)ft(xf)h = \{n + l)h n (x n )5(x - x')} , 

ffl(x')t/j(x)h = nh n (x n ~ 1 , x)S(x n — x') 

Suppose that we are considering a subspace of (1.51) whose vectors consist of elements 
having some fixed symmetry, for instance, the functions h n (x n ) are symmetric or anti- 
symmetric. But, as is seen from above, the action of field operators does not preserve 
this symmetry, except for the action of the annihilation operators, ifj(x)ip(x'). This 
happens because the domain of the field operators, when this domain is a subspace of 
(1.51), generally, does not coincide with their range. Recall that by the range of field 
operators we mean the range of operators (1.56). 

It is clear, that to obtain simple commutation relations, one needs to invoke some 
additional conditions. Before doing this, let us explain why (1.54) and (1.55) are called 
annihilation and creation operators, respectively. 

Choose as a basis of (1.51) that consisting of vectors 








\n) := [S mn h n (x n )] = 





K{x n ) 






(1.57) 



A vector adjoint to (1.57) is a row 

(n\ := [0 ... h* n (x n ) 00 . . .] . 

The notation used in (1.57) and (1.58) is called the Dime notation. 
The scalar product of elements from the basis {\n)} is 

(m\n) = S mn (h n , K). 
So, it is a complete orthogonal basis. It is easily normalizable by putting 

(n\n) = (h n ,h n ) = I. 
Then, (1.57) is called the n-particle state. 



1.58) 
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The O-particle state, that is 



|0) := [5 m0 h } 



h 





1.59) 



with \h \ = 1, is called the vacuum vector or vacuum state. 

The action of the annihilation operator (1.54) on the vacuum vector (1.59) is 

$(x)\0)=0, (1.60) 

and its action on other vectors from the basis with n — 1, 2, . . ., gives 



tp(x)\n) = ^fn 5 mn _ih n {x n 1 ,x) 



(1.61) 



The action of the creation operator (1.55) on an n-particle vector (1.57) yields 



ft(x)\n) = Vn + 1 [S mn+1 h n (x n )5{x n+1 - x)} . 
For the operators (1.56) we have from (1.60) and (1.61) 

1?/|0> = o, 



(1.62) 



with 



ti n „ 1 {x n ~ 1 )= J h n (x n -\x)f(x)dx, 



which shows that ipf transforms an n-particle vector to an (n — l)-particle vector. In 
other words, t/jf annihilates a particle. From (1.62) we get 



with 



ty(x)\n) = VnTT S mn+1 h' n+1 (x 



h'(x n+1 )=h n (x n )f(x n+1 ), 



which shows that tp^ transforms an n-particle vector to an (n+ l)-particle vector, that 
is, creates a particle. 

In the same way one can check that the double action of the field operators, given 

by 

ijj(x)ft(x')\n) = (n + 1) [S mn h n (x n )S(x - x')} , 
tp\x')ip(x)\n) = n 5 mn h n (x n ~ 1 ,x)S(x n - x') , 
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means that starting from an n-particle vector and first creating a particle and then 
annihilating it, or first annihilating a particle and then creating it, one always returns 
to an n-particle vector, though not necessarily the same. 

By acting n times on the vacuum state (1.59), with h = 1, by the creation operator, 

n 

i=l L i=l 



II^Vi)|0> = Vn! 



one can construct any n-particle vector (1.57) as 

\n) = ^=l h n (x n )f[ft( Xi )\0)dx n (1.63) 
and any vector of space (1.51) as 

oo 

h=J2 C n \n), 



n=0 



where C n are complex numbers. 

Till now, we have not required from the vectors of the space (1.51) any symmetry 
properties. But, according to the symmetry principle discussed in Section 1.1, the 
quantum states of identical particles are to be either symmetric or antisymmetric with 
respect to permutation of their variables. Symmetric quantum states correspond to 
bosons, antisymmetric to fermions. Thus, we have to construct the subspace of (1.51) 
possessing the desired symmetry properties. 

The set x n of the variables xi, x 2 , ■ ■ ■ , x n permits one to make n! permutations. 
Enumerate these permutations by the index % = 0, 1, 2, . . . , (n — 1)! and denote the 
corresponding permutation operators by Pi(x n ), so that 

P (x n ) = 1, P?(x n ) = 1. 

Using arbitrary functions h n (x n ), we may construct symmetric and antisymmetric 
functions by the receipt 

= ^'(^TPiixnhnixn. (1-64) 
Vn! i=0 

By construction, from (1.64) it follows 

Pi{x n )f±{x n ) = ±f±(x n ) (i > 1). (1.65) 
Denote a Hilbert space of such functions as 

K = {f±{x n )}. 

Analogously to (1.51), define the Hilbert space 

n ± := ©~ H± (1.66) 
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of vector columns 

f± = ifH* n )} e n ± . (i.67) 

The symmetric space 7i + represents Bose states while the antisymmetric space 7Y_, 
Fermi states. A space having the structure (1.66), possessing a vacuum vector, is a 
Fock space. 

The symmetric and antisymmetric spaces in (1.66) are subspaces of (1.51), and can 
be considered as projections 

H± = P±H (1.68) 

of H by means of projectors P + and P_. 
Introduce the field operators 

ipi(x) := P±ip t (x)P±, (1.69) 

such that ipf + and i>f + , defined as in (1.56), always act inside the symmetric space 
H+, while ipf- and i/)f_, inside the antisymmetric space 7i-. By this definition 

tp+(x)^(x r ) = ip + (x)ipl(a/) = 0. (1.70) 

In what follows, for simplicity of notation, we shall write ip(x), representing either 
ip+(x) or if)-(x) from (1.69), and ip\x) being either ip + (x) or ip^(x). Respectively, 
we will write / instead of f± from (1.67). To distinguish a commutator from an 
anticommutator, we denote them as 

[A,B] T := ABtBA. 

Using the formulas of this Section for the field operator (1.54) and (1.55), for the 
operators (1.69) we obtain the commutation relations 

[^(x),^(x')] T / = 0, 

ty(x),^(x')] T f = 5(x-x>)f, (1.71) 

the upper sign being related to Bose- and the lower, to Fermi-cases. One also says 
that the commutators in (1.71) define Bose statistics while the anticommutator, Fermi 
statistics. Since the relations (1.71) are valid for all corresponding vectors from (1.66), 
a vector / in (1.71) can be omitted, for example 

[ip(x),^(x')} T = 5(x-x')l, 

where the unity matrix is 

i:= [5 mn 5{x n -y n )\. 
Usually, for brevity, one omits the unity matrix as well. 
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Recall again that ip(x) and ij)\x) are operator distributions. Actual field operators 
are given, as in (1.56) by 

*p f = J f(x)ip(x)dx, ip} = J f(x)^(x)dx. (1.72) 

The commutation relations for (1.72) are 

[l^]=F = [l^J]=F = o, 

for any square-integrable functions f(x) and g(x). 

The time evolution of the field operators is governed by a Hamiltonian H which 
is a matrix of the form (1.53). The evolution operator U(t) is defined, analogously to 
(1.43) and (1.44), as a unitary operator satisfying the Schrodinger equation 

ih^U{t) = HU(t) (1.73) 

and possessing the properties 

U(t)U(t') = U(t + t'), 

U+(t) = U- 1 (t) = U(-t), (1.74) 
U(0) = I. 

An annihilation field operator at time t is given, as in (1.46), by the relation 

ijj(x,t) := U + (t)iP(x,0)U(t), (1.75) 

in which 

ij}(x,0) :=il>(x). (1.76) 

In the same way, the creation field operator ^(x,t) is defined. 
The evolution equation for ip(x,t) follows from (1.75), 

ih^(x,t) = [il>(x,t),H(t)] (1.77) 

where [, ] means a commutator and 

H(t) := U + {t)HU{t). 

The equation (1.77), together with its Hermitian conjugate giving the evolution equa- 
tion for ip^(x,t), are the Heisenberg equations for the field operators. 
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1.4 Local Observables 



The operators corresponding to observable quantities are constructed with the help of 
the field operators. In doing this, one has to keep in mind that the manifold V can 
be unbounded. Therefore, first one has to consider submanifolds of V which are open 
bounded manifolds. Let A C V be such a submanifold. The general structure of an 
operator of a local observable is 

00 1 r 

A A := ]T -7==/ ^(x m )...^(x 1 )A(x m , y n )Hyi) ■ ■ ■ i>(y n )dx m dy n , (1.78) 
m ,n=o Vmlnl J A 

where A(x m ,y n ) can be a generalized function or an operator defined on functions 
with the domain A m x A n . Integration in (1.78) for all variables goes over the same 
submanifold A. 

For example, the local number-of-particles operator is 

N A = [ ^(x)ip(x)dx. (1.79) 

J A 

Using the density-of-particle operator 

n(x) := i>\x)i>(x), (1.80) 

(1.79) can be written as 

N A = J h{x)dx. 
The action of (1.80) on an n-particle vector (1.57) gives 

n(x)\n) = n [5 mn f n {x n ~ 1 , x)5{x n - x)\ , 

from where 

N\\n) = n\n). 

The n-particle vectors are eigenvectors of N A , the numbers n — 0,1,2, . . ., being eigen- 
values. This explains why (1.79) is called the number-of-particle operator. 
The time evolution of operators (1.78) follows from (1.75) giving 

A A (t) = u + (t)A A u(t). (i.8i) 

The evolution equation for (1.81), in compliance with (1.73), is the Heisenberg equation 

ihj t Mt) = i^A A (t) + [A A (t),H(t)\, (1.82) 



where 



a a a 

-A A (t):=U + (t)^U(t), (1.83) 



which is equivalent to (1.47) and (1.48). 
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Introduce the statistical operator, or density operator, p(t), which is a positive 
bounded self-adjoint operator, normalized to unity, 

Tr p(t) = 1, (1.84) 

and satisfying the Liouville equation 

ihj t p{t) = [H{t),p(t)]- (1-85) 

The time evolution, following from (1.85), is 

P {t) = U(t)pU + (t); p:=p(0). (1.86) 

The mathematical expectations for the operators of local observables are called the 
statistical averages and are defined as 

(i A > := Tr p{t)A A = Tr pA A (t), (1.87) 

which is analogous to (1.50). The trace in (1.87), as well as in (1.84), is taken over any 
basis from either 7i + or according to the statistics of particles. Since the average 
(1.87) depends on time, we may write 

(Aa) = (A A (t)). 

Define the statistical density matrices through their elements 

p mn (x m , y n , t) = Tr ^(xi) . . . ^(x m )p(t)^(y n ) . . . ^( Vl ). (1.88) 

Then the average (1.87) takes the form 

00 1 r 

(A A )= ]T -7==/ A{x m ,y n )p mn (y n ,x m ,t)dx m dy n . (1.89) 
m , n =o Vmlnl J A 

The widespread situation is when 

A(x m ,y n ) = 5 mn A(x n ,y n ). 

Then we need only the elements (1.88) that are diagonal with respect to the indices 
m, n, that is 

Pn {x n ,y n ,t):=p rm {x n ,y n ,t). (1.90) 

In this case, (1.89) reads 



n=0 

A linear space 



00 1 r 

<^a> = E -T / A ( xn , V n )Pn(y n , x n , t)dx n , dy n . (1.91) 
A A := {A A \ A C V} (1.92) 
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of the self-adjoint operators of local observables is an algebra of local observables [26,27]. 
A set 

(A A ) := {<i A », (1.93) 

formed by the statistical averages (1.89), is called the statistical state. 

If the manifold V is unbounded, we can define an ordered sequence {Ai C A 2 C 
. . . C V} of open submanifolds Ak, such that 

lim A k = V, 

k— +oo 

and define a set of algebras, {A Al C A A , 2 C . . . A Ak C . . .}. The inductive limit 

Ay := lim A Ak (1.94) 

k—>oo 

is called the quasilocal algebra, or global algebra of observables. 

1.5 Order Indices 

The statistical density matrices 

p n (t):=[p n (x n ,y n ,t)} ((*>:= 1) (1.95) 
with elements (1.90) are matrices in the x — y space. Their action on a vector 

fn := [f(x n )], 

being a column with respect to x, is given by the rule 



Pn{t)fn ■-- 



p n (x n ,y n ,t)f n (y n W 



The elements (1.90), owing to the evolution properties (1.75) and (1.86), can be written 
in two forms: 

p n (x n , y n , t) := Tr ^{x 1 ) . . . ^{x n )p{t)^{y n ) . . . ^( Vl ) = 

= Tr i/>( Xl , t)... ij>(x n , t)p^{y n , t)... ^( Vl , t). (1.96) 
According to the commutation relations (1.71), the symmetry property 

P(x i ,x J )p n (x n ,y n ,t) = ± Pn (x n ,y n ,t) 

holds. 

The elements (1.96) cannot, in general, be factorized into products of functions 
depending separately on x n and y n . This means that, by construction, the statistical 
state is a mixed state. For latter, as is discussed in Section 1.2, all density matrices 
with n — 0, 1, 2, . . . can be called reduced. This is especially natural when working in 
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the field-operator representation where the number of particles is not fixed. Then one 
calls p n (t) simply an n-matrix. 

However, we may fix the average number of particles 

N := (N v ) = J Pl (x,x,t)dx, (1.97) 

where the integration is over V, and 

x, t) = (n(x)) = Tr p(t)n(x). (1.98) 

Then, it would be reasonable to distinguish among all n-matrices those for which 
n < N and to apply the name reduced only to them. 

Let us fix the average number of particles (1.97) and consider a subspace of the 
Fock space for which 

(N V A V ) = N(A V ). (1.99) 
Then for the 2-matrix with elements 

p 2 (x 1 ,x 2 , 2/1,2/2, t) = {il)\y2)^{yM{xM{x2)), (1.100) 

using the commutation relations (1.71), we can write 

P2(xi,x 2 ,yi,y 2 ,t) = (^(y 2 )ip(x 2 )ip 1 (yi)ip(x l )} - 5(x 2 - yi)(^ t (y 2 )V'(^i)>- 
From here, taking account of (1.99), we get 

J p 2 (xi,x, x 2 ,x, t)dx — (N — l)pi(xi,x 2 ,t). (1.101) 

Integrating (1.101) once more, and using (1.97), we find 

[ p 2 (x 1 ,x 2 ,x 1 ,x 2 ,t)dx 1 dx 2 = N(N - 1). 

Generally, the following recursion relation takes place: 

J p n (x n -\x, y n ~\ x, t)dx =(N-n+ l)p n ^{x n ~\y n -\ t). (1.102) 

The normalization condition for an n-matrix is 

J p n (x n ,x n ,t)dx n = {N m n)l ■ (1-103) 

Note that the normalization (1.103) for p n (t) differs from the normalization (1.27) 
for D n (t). 

Consider a system with a fixed average number of particles N and occupying in the 
real space a volume V. One says that the thermodynamic limit is performed if 

N 

iV->oo, V ^ oo , — -> const , (1.104) 
24 



which is denoted as limjv^oo or lirriy-^oo- 

An ensemble of particles is called a thermodynamic system, if the statistical state 
(1.93) under the thermodynamic limit (1.104) behaves so that 



lim 

N^oo 



< oo. (1.105) 



For example, the limit 

p:= lim Ln v ) (1.106) 
k— >oc y 

means the average density of particles, and the limit 

lim ±-E N := lim 1<£ V ) (1.107) 

gives the internal energy per particle. 

To understand the behaviour of density matrices under the thermodynamic limit, 
define the norm 

||p„(*)||:= sup ( Pn (t)f njPn (t)f n ) l/ \ (1.108) 

ll/n||<l 

where p n {t) is the density matrix (1.95) and 

| \fn\ | - = {fn, fn) ^ ■ 

The statistical operator p(t) is, by definition, a positive self-adjoint operator. Then, the 
density matrix (1.95) with elements (1.96) is also positive and self-adjoint. Therefore, 
the norm (1.108) can also be presented as 

\\p n (t)\ \ = SUp \(f n ,Pn(t)fn)\ = 

\\fn\\<l 

= sup / f*(x n )p n (x n ,y n ,t)f n (y n )dx n dy n . (1.109) 
Define the order indices [28-31] for density matrices by the limit 

^w-J&nl^r 1 - (1110) 

where the thermodynamic limit (1.104) is assumed and n — 0, 1,2, . . .. Each density 
matrix p n (t) is positive and obeys the normalization condition (1.103), because of 
which, using the Stirling formula for iV — > oo, we have 

0<a n <n. (1-111) 

If the eigenvalue problem 

Pn(t)(f nq (t) = 7 nq (t)(p nq (t) 
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is solvable for an n-matrix, then 



||p„(t)|| = sup7„ g (t) := 7 n (i), 
<? 

and the order indices (1.110) can be expressed through the largest eigenvalues of the 
corresponding n-matrices: 

This shows the importance of the largest eigenvalues of density matrices in the 
thermodynamic limit. 

1.6 Information Functionals 

To define the density matrices with elements (1.96), we need to know the statistical 
operator pit). The latter obeys the Liouville equation (1.85). This operator equation 
not always defines p(t) uniquely but some additional conditions may be required. Such 
conditions can be formulated by invoking the notions of information theory [32]. 

In what follows, implying that the operators of local observables are associated 
with a manifold V, for simplicity, we shall not ascribe the corresponding index at A v 
writing instead just A. Also, remind that if A is an operator from the algebra A of 
local observables, and FL4(i)] is a function of A(t), then the evolution of FL4(£)] is 
given by the law 

F[A(t)\ = F [u + (t)Au(t)} = U + (t)F[A]U(t), 

where A := A(0). This law can be considered either as a definition or it can be proved 
if the function F[A] is defined as a Taylor expansion in powers of A. Similarly, if F[p(t)} 
is a function of the statistical operator, then, in compliance with (1.86), 

F[p(t)} = F [u(t)pU + (t)\ = U(t)F[p]U + (t), 

where p = p(0). 

Introduce the information operator 

I{t) := \np(t). (1.113) 

The time evolution for (1.113), owing to the evolution laws discussed above, is 

i(t) = u(t)i(o)u + (t), (i.ii4) 

where 

1(0) := lnp(0) = In p. 
Consequently, the information operator (1.113) satisfies the Liouville equation 

ihj t i(t) = [H(t),i(t)]. (i.ii5) 
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The Shannon information is defined as a functional 

I s (t) := Tr p(t)I(t) = Tr p(t) lnp(t). (1.116) 

Its relation with 

Is ■= ls(0) =Trp In p (1.117) 

easily follows from the equality 

Tr p(t)I(t) = Tr p(0)/(0) 

resulting from (1.86) and (1.114). Therefore, the Shannon information (1.116), actually, 
does not depend on time, 

I s (t) = I s (0). (1.118) 

The equality (1.118) is the essence of the Liouville theorem. 
The information entropy is given by 

S := -I s = -Tr phi p. (1.119) 

Since the statistical operator p is a positive bounded operator normalized to unity by 
(1.84), the information entropy (1.119) is positive and, because of (1.118), it has no 
time evolution. 

It is important to emphasize that the Shannon information (1.116) is not an average 
of the information operator (1.113). An average information can be defined as 

I A (t) := Tr pl(t) = Tr p\np(t). (1.120) 

Its time evolution follows from the Liouville equation (1.115) giving 



ihj t I A {t) = Tr p[H{t)J{t) 



[1.121] 



while the Shannon information, in comformity with (1.118), is constant in time. One 
may expect that the average information (1.120) decreases with time because of the 
following. 

Proposition 1.1. (Nonincrease of average information). 

Let the average information be given by the functional (1.120), and let an inverse 
operator p _1 (t) exist satisfying 

p(t)P~ 1 (t) = l, 

then for t > 

I A (t) < I A (0). (1.122) 

Proof: Consider the difference 

I A (0) - I A (t) = Trp [lnp-lnp(t)]. 
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Due to the existence of the inverse operator f>~ l {t), we can write 

In p — In p{t) = In pp^ 1 (t) . 

Since, by definition, p(t) is a positive bounded operator in a Hilbert space, its inverse 
is a positive operator. The product of two positive operators is a positive operator, 
so pp~ l {t) is a positive operator. For a positive operator I in a Hilbert space the 
inequality 

lnX^l-X" 1 (X>0) 

holds, if the inverse operator exists. Applying this inequality to the positive operator 
pp~ l {t), we have 

\npp-\t)>l-p' l p{t). 

Consequently, 

I A (0) - I A (t) >Tr p[l- p- l p(t) . 

As far as the statistical operator p(t) for any t > is normalized to unity by condition 
(1.84), we get 

Tr p 1 — p~ l p{t) = Tr p — Trp(t) = 0. 
Thus, we find that Tr pin p > Trp In pit) and 

I A (0) - I A (t) > 0, 

which prooves (1.122). 

If, instead of the entropy (1.119), we introduce the average entropy 

S A (t):=-Trp\np(t), (1.123) 

then from the above proposition, it follows 

S A (t) > S A (0) it > 0). (1.124) 

The Shannon information (1.116) is an unconditional information. It may happen 
that some additional knowledge about the properties of a considered system is available. 
Then conditional information functionals can be defined. 

Assume that for a set {Aj(t)| i — 0, 1, 2, . . . , n} of operators from the algebra A of 
local observables some statistical averages 

A^t) := (Ai(t)) = Tr pA.it) = Tr p(t)Ai (1.125) 

are known. One obvious case is for the unity operator A^t) — 1, for which we know 
that 

(i) = Tr pit) = 1. (1.126) 
Introduce the operator of conditional information 

n 

I con it) :=]np(t) + ^\i(t) [A t -A^t)] (1.127) 

i=0 
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and the conditional information 

Iconit) := Tr p(t)I con (t) = 

n 

= I s + £ K(t) [Tr pit) A, - Ai(t)] , (1.128) 

i=0 

where Xt(t) are Lagrange multipliers satisfying the variational condition 

SIcon(t) 



SXi(t) 



0. (1.129) 



The meaning of (1.129) is to guarantee the validity of (1.125). Note that, contrary to 
the Shannon information (1.116) that does not depend on time, according to (1.118), 
the conditional information (1.128), generally, is a function of time. 

If the behaviour of some observable quantities, characterized by a set {Ai(x, t)\ i — 
0, 1, 2, . . . , n} of operators Ai(x, t) is assumed to be known in each local vicinity of x 
and t, this is described by the averages 

Ai(x, t) := (Mx, t))=Tr p{t)A\{x). (1.130) 

In such a case, one can define the operator of local conditional information 

n „ 

I loc (t) := lnp(t) / \(x,t) [Mx) - A{x,t)] dx (1.131) 

i=o J 

and, respectively, the local conditional information 

I loc (t) :=Tr p(t)I loc (t), (1.132) 
in which the variational condition 

SIloc(t) 



SXi(x,t) 



= (1.133) 



defines the Lagrange multipliers Xi(x,t) guaranteeing the validity of (1.130). The local 
conditional information (1.132) also depends on time. 
Finally, when the transition amplitudes 

Ai(x, t, x', t') := (Ai(x, t, x', 0) (1-134) 

are given as the averages of operators A^x, t, x', t') with the evolution property 

Ai{x,t,x',t') = U + (t)Ai(x,0,x',t')U(t), 

then one can introduce the operator of nonlocal conditional information 



n . 

I non (t) := lnp(t) +]T / Xi{x,t,x',t') 

i=0 



X 
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x Ai(x,0,x',t') - Ai(x,t,x',t') dxdx'dt' (1.135) 
and the nonlocal conditional information 

I non (t) := Tr p(t)I non (t) (1.136) 

with the Lagrange multipliers defined as in (1.129) and (1.133). Undoubtedly, (1.136) 
is a function of time. 

When the Liouville equation (1.85) does not allow a unique definition of the statisti- 
cal operator, then one can invoke the principle of minimal information stating that the 
statistical operator must have such a form that provides the minimum of information. 
This implies the validity of the extremum condition 

in which I(t) is an information functional for either the conditional information (1.128), 
local conditional information (1.132), or nonlocal conditional information (1.136), de- 
pending on a case considered. Also, the inequality 

Si- 

is to be true, meaning that the left-hand-side of (1.138) is a positive operator. The 
corresponding Lagrange multipliers are given by the variational conditions as in (1.129) 
or (1.133). Performing the variation of information with respect to the statistical 
operator p(t), one assumes the existence of the inverse operator p _1 (t), so that 

\n P {t)=p-\t), 



5p(t) 



and one uses the equality 



5p(t) 



Tr p(t)A(t) = A(t). 



For any kind of a statistical operator, the normalization condition (1.126) is sup- 
posed to hold. So, one of the Lagrange multipliers, namely A (t), should have the 
same form for all statistical operators. This Lagrange multiplier may be related to the 
partition function Z{t) through the equation 

Z(t) = exp{l + A (t)} . (1.139) 

In the case of the conditional information (1.128), the principle (1.137) yields 



Pit) = {- E (1-140) 
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where the partition function 

Z(t)=Trexp\-f2WtM (1-141) 



i=i 



plays the role of a normalization factor guaranteeing the validity of (1.126). The 
inequality (1.138) holds true since it tells that 

iM =p (t)>0 ' 

which is actually the far as p(t) is a positive operator. 

The local conditional information (1.132), in compliance with the extremum equa- 
tion (1.137), leads to the local statistical operator 

Piodt) = ^^yexp j-]T J X l (x,t)Mx)dx^ (1.142) 

with the local partition function 

Z loc (t) = Tr exp j-fj/ \(x,t)A i (x)dx^ . (1.143) 

In the same way, for the nonlocal conditional information (1.136) one finds the 
nonlocal statistical operator 

Pnon{t) = - ~~7~t exp ( - [ X i( x i *, ^ 04(s, 0, x' , t')dxdx ' dt' 1 (1.144) 

^non{t) { i= l J J 

with the nonlocal partition function 

Znonit) =Tr exp j-^y \ i (x,ty,lf)A i (x,Q,a?,t')dxdx , dlf'^. (1.145) 

Although the found statistical operators p(t), pi oc (t), and p ncm {t) depend on time, 
however this dependence does not necessarily coincide with the time evolution pre- 
scribed by the Liouville equation (1.85). To correct this deficiency, one can do the 
following. Introduce the Liouville super-operator C, acting on the statistical operator 
pit) by the rule 

tp{t)--\[H{t),p{t)], (1.146) 

and the operator 

6:=(±-it) m . (1.147) 
Define the operator of generalized conditional information 

Igen(t) 1= Ln(t) + \(t)6 (1.148) 
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and the generalized conditional information 

Igenit) 1= I con (t) + \{t)Tr f>(t)6. 



1.149) 



Using the principle of minimal information, we obtain the generalized statistical oper- 
ator 

= v^m ex P |- E Ht)A - \(t)d) (i.i5o) 

with the generalized partition function 



Z gen (t) = Tr exp - £ A,(t)i, - \{t)6 



(1.151) 



i=i 



The Lagrange multiplier X(i) is to be defined from the variational condition 

SIgen(t) 



5\{t) 



0, 



which yields 



Trp(t)6 = Trp(t) (±-iAp(t) = 0. 



The latter equation, together with (1.146) and (1.147), gives 



Tr pit) I ^i- pit) -[H{t), pit) 



0, 



(1.152) 



that is the averaged Liouville equation (1.85). In the same way, one may correct the 
local statistical operator (1.142) and nonlocal statistical operator (1.144). 



1.7 Equilibrium States 

The statistical state (A), which is the set {(A(t))} of expectation values (A(t)) of the 
operators A(t) representing observables quantities, generally, depends on time. It may 
happen, that this dependence is absent, so that 

(A(t)) = (A(0)) . (1.153) 

Then (1.153) presents an equilibrium state, or time-invariant state. 

Since (1.153) must hold for all operators from the algebra of observables A, the 
necessary and sufficient condition for (1.153) is the time invariance of the equilibrium 
statistical operator, 

pit) = p(0) := p. (1.154) 

This means that 

j t p{t) = 0, (1.155) 
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and the Liouville equation (1.85) yields 

if(t),p(i)] = 0. (1.156) 

Eq.(1.156) is valid if and only if the statistical operator in (1.154) is a function of H{t) 
and of operators Ci(t), with i — 1, 2, . . . , n, for which 

Ci(t),H(t)]=0, J^i(*) = 0. (1.157) 

The operators satisfying (1.157) are called the integrals of motion, because, according 
to (1.82) and (1.157), 

jC.it) = 0, Ci(t) = (7,(0) := Ci . 

From (1.154) and (1.156) it is evident that the Hamiltonian 

H(t) = H(0) := H := C 1 (1.158) 

is one of such integrals of motion. The expectation value 

E := (H) (1.159) 

is called the internal energy of the system. 

The knowledge that the equilibrium statistical operator (1.154) is a function of 
integrals of motion is not sufficient for defining its concrete form. The latter can 
be found with the help of the principle of minimal information as described in the 
previous Section. To this end, one constructs the conditional information (1.128) with 
the operators A. identified with integrals of motion, 

Ai = Ci (i = l,2,...,n); A = i . 

As a result, one comes to (1.140) and (1.141) but without time dependence, 

p=lexp(-X>A), (1-160) 

where 

Z = Tr exp^-X>C^ . (1.161) 

The statistical operator of the form (1.160) is named the Gibbs statistical operator 
and the statistical state (A), composed of averages (A) defined as Tr pA with this 
statistical operator, is a Gibbs state [33-35]. 

There can be different Gibbs states depending on the given integrals of motion. 
For example, if the sole integral of motion is a Hamiltonian H, then the conditional 
information (1.128) becomes 

I con = Tr plnp + A (Tr p - 1) + Ai (Tr pH - E) , 
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where 

A :=lnZ-l, Ai := /? . 
The principal of minimal information gives 

p = le-^ (1.162) 

with 

Z = Tr e~ pA . (1.163) 

The parameter (5 is called the inverse temperature, and it is related to the absolute 
temperature measured either in energy units, 6, or in Kelvin degrees, T, by the equa- 
tions 

pe = i, e = k B r, 

in which k B is the Boltzmann constant. 

The statistical operator (1.162) is named the Gibbs canonical operator. One says 
that a physical system which it describes is represented by the Gibbs canonical ensem- 
ble. 

A quantity F satisfying the equation 

Tr exp (- £ = exp {-(3F) (1.164) 

is termed the thermodynamic potential. Another form of (1.164) is the definition 

F : = -0 InTr exp ^- J2 A ^ ■ (1.165) 

The thermodynamic potential for the Gibbs canonical ensemble, according to 
(1.163) and (1.165), is 

F := -0 In Z = -0 In Tr e-P 6 . (1.166) 

For a system of N particles inside a real-space volume V, the thermodynamic potential 
(1.166) which is a function 

F = F(Q, V, N) (1.167) 

of temperature, volume, and number of particles, is called the free energy. 
Writing for (1.167) the differential form 

dF = -S c dB -p c dV + fi c dN, (1.168) 

we define the canonical entropy 
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canonical pressure 

Pc ■= - (%) , (1-170) 



and the canonical chemical potential 



H3e,' (Lm) 

Consider now the case of two integrals of motion, one of which is the Hamiltonian 
H and the second is the number-of-particles operator 

N — N(t); [N,H]=0. (1.172) 

The average number of particles is 

N:=(N)=TrpN. (1.173) 

As the conditional information (1.128) we have 

I con = Tr plnp + A (Tr p - 1) + A x (TV pH - E) + A 2 (Tr pN - n) , 

in which 

A :=mZ-l, X 1 :—P, A 2 := -fin • 
The principle of minimal information yields 

p =L e -m-^) (i.i74) 
z 

with the partition function 

Z = Tr e-ft*-^ . (1.175) 

The parameter p introduced above is the chemical potential. 

The statistical operator (1.174) is called the grand canonical operator, and one says 
that a physical system described by it is represented by the grand canonical ensemble. 

For the thermodynamic potential (1.165), we get 

D, := -0 In Tr e -^ A ~^ , (1.176) 
which is termed the grand potential. The latter is a function 

n = n(e,v,fi) (1.177) 

of temperature, volume, and chemical potential. 
The differential form 

dVL = -S g d& - pgdV - N g dp ( 1 . 1 78 ) 
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for (1.177) defines the grand entropy 




dQ 



) 



(1.179) 



the grand pressure 




) 



(1.180) 



and the grand number of particles 




dtt 



) 



e,v 



(1.181) 



Similarly, one may construct other statistical operators and the corresponding thermo- 
dynamic potentials. 

The quantities obtained by differentiating thermodynamic potentials are named 
thermodynamic characteristics, or thermodynamic functions, and their variables, ther- 
modynamic parameters, or thermodynamic variables. 

Different thermodynamic characteristics defined for different ensembles are, of 
course, different functions not necessary having simple relations with each other. This 
is because, generally, different ensembles represent different physical systems. An en- 
semble correctly representing a given physical system is a representative ensemble. 

It may happen that the same system of particles is equally representable by several 
ensembles, which means that the analogous thermodynamic characteristics of different 
ensembles can be transformed to each other by a corresponding change of thermody- 
namic variables. 

Representative ensembles equally describing the thermodynamic properties of a 
physical system are thermodynamically equivalent ensembles. 

For example, if the canonical and grand canonical ensembles are thermodynamically 
equivalent, then the canonical entropy (1.169) and the grand entropy (1.179) are related 
with each other by the equations 



the canonical pressure (1.170) and the grand pressure (1.180) satisfy the relations 



and the relations between the canonical chemical potential (1.171) and the chemical 
potential //, and between the grand number of particles (1.181) and the number of 
particles N are 



s c (e,v,N g (e,v,fi)) = s g (e,v,fx), 



S g (Q, V, // c (9, V, N)) — S C (G, V, N) 



Pc (e,v,N g (e,v,ri)= Pg (e,v,fx), 



Pg (e,v,fi c (e,v,N)) = Pc (Q,v,N) 



fx c (Q,V,N g (e,V, = 
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N g (Q,V,t, c (e,V,N)) = N . 

In such a case the differential forms (1.168) and (1.178) are related with each other by 
a Legendre transformation which gives 

n(e, v, /i) = f(6, v, N g (e, v, /i)) - fiN g (e, v, //), 

F(0, V, iV) = fi(6, V, /i c (0, V, iV)) + ^(6, V, N)N. 
Keeping this change of variables in mind, it is customary to write the latter relation as 

Q = F-fiN. (1.182) 

Traditionally, one also omits the indices at the thermodynamic characteristics writing 



V,N 





"=l^» N = 



(1.183) 

9, ft 

e,v \^/ey 
and employing in each particular case the corresponding dependence on thermodynamic 
variables. 

For real physical systems the thermodynamic equivalence of ensembles usually oc- 
curs only in the thermodynamic limit (1.104). Then the above relations involving the 
change of thermodynamic variables are to be considered as asymptotic equalities. 

The thermodynamic equivalence of ensembles does not necessarily involve the equiv- 
alence of Gibbs states. This means the following. Suppose the canonical ensemble is 
thermodynamically equivalent to the grand canonical ensemble resulting in (1.182) and 
(1.183). However, the set of averages for the operators from the algebra of observables, 
calculated with the Gibbs canonical operator (1.162) does not necessarily coincide with 
the set calculated with the grand canonical operator (1.174). 

Let us denote by (A) c and (A) g the Gibbs states defined through the Gibbs canon- 
ical operator (1.162) and the grand canonical operator (1.174), respectively. We shall 
write 

(A) c =(A) g , (1.184) 

if all corresponding averages coincide with each other. Call (A) c the Gibbs canoni- 
cal state and (A) g the grand canonical state. Then the equivalence of the statistical 
states (1.184) is termed the statistical equivalence of ensembles. The thermodynamic 
equivalence follows from the latter, but not reverse. 

The definition of entropy given for the Gibbs canonical state in (1.169) or for the 
grand canonical state in (1.179) can be generalized for an arbitrary Gibbs state with a 
thermodynamic potential (1.165) as 

dF 

S:=-qq, (1-185) 
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where the partial derivative with respect to temperature is taken under all other ther- 
modynamic variables kept fixed. So defined entropy (1.185) is the thermodynamic en- 
tropy. Its particular examples are the canonical entropy (1.169) and the grand entropy 
(1.170). We may remember also that we have had the definition of the information 
entropy (1.119). What are conditions under which the thermodynamic and information 
entropies coincide? 

Proposition 1.2. (Criterion for entropy equivalence). 

The thermodynamic entropy (1.185) for a Gibbs state with the Gibbs statistical 
operator (1.160) coincides with the information entropy (1.119), that is 

OF 

- — = -Trplnp, (1.186) 
if and only if the Lagrange multipliers satisfy the equality 

g (§H) ^ = ■ (L187) 

Proof: From the definition of the thermodynamic entropy (1.185), with the ther- 
modynamic potential (1.165), we have 

S = ^{t^(C i )-Fj . (1.188) 

Now, let us consider the information entropy (1.119). For the statistical operator 
(1.160), we get 

n 

lnp = -^AiCi-lnZ . (1.189) 
i=i 

With the thermodynamic potential F defined in (1.166), this results in 

- Tr plnp = (3 ^(C t ) - F^j . (1.190) 

The comparison of expressions (1.188) and (1.190) shows that they coincide if and only 
if the condition (1.187) holds, which concludes the proof. 

Corollary. The canonical entropy (1.169) and the grand entropy (1.179) coincide 
with the information entropy (1.119) expressed in the corresponding variables. 

In the case of the Gibbs canonical state, 

and the criterion (1.187) is valid. The thermodynamic entropy has the form 

S — (5(E — F) 



38 



coinciding with the information entropy. 

For the grand canonical state the Lagrange multipliers are 

Ai = (3, \ 2 = -fin 
and the criterion (1.187) holds since 

<9Ai _ 9A2 

The thermodynamic entropy (1.188) equals the information entropy (1.119) being 

S = [3(E — fiN — Q) . 
Note that in these two cases we have the equality 

^1 = ^1 (1191} 

d/3 (3 ' 1 > 

which is a sufficient condition for (1.187) 

Similarly to the considered cases, we may define other statistical states described 
by the Gibbs statistical operator (1.160) with other integrals of motion. For example, if 
the real-space volume of the system is not fixed, we can introduce the volume operator 
V, which is a projector from the real space to a chosen volume, such that 

(V) = V. (1.192) 

Then, constructing the corresponding conditional information and minimizing it, one 

comes to the Gibbs isobaric state with the Gibbs isobaric operator 

p=±e-n*-M, (1.193) 

Zj 

in which 

Z = Tr e-^-v*) . 
The thermodynamic potential related to the Gibbs isobaric state, 

G := -6 InTr e'^ 6 ^ , (1.194) 

is called the Gibbs potential. It is a function 

G = G(0,p,N) (1.195) 

of temperature, pressure and number of particles, with the differential form 

dG = -SdQ + Vdp + ndN . (1.196) 
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The later defines the isobaric entropy 



the average volume 

'dG\ 



V -[^)e,' <L198) 



and the isobaric chemical potential 

,==(§£) • d.l») 

\ / e,p 

In formulas (1.196) to (1.199), following tradition, the index specifying the affiliation of 
the thermodynamic characteristics to an ensemble, in this case to the isobaric ensemble, 
for brevity, is omitted. The relation of the Gibbs potential (1.194) with the free energy 
and the grand potential is 

G = F + pV = Vt + /jlN + pV . (1.200) 

The general form of the thermodynamic potential (1.165) has remarkable convexity 
properties with respect to Lagrange multipliers. The first derivative of (1.165) is 

dF 



Differentiating this once more, we find 

d 2 F 



-e (<c?> - (6^) 



0A? 

For any self-adjoint operator A, the quantity 

(A 2 ) - (i) 2 = <(i - (A)) 2 ) 



is nonnegative. If integrals of motion pertain to the algebra of observables, they are 
self-adjoint. In such a case, 

d 2 F 

^<0 (i>l), (1.201) 

that is, the thermodynamic potential is a concave function with respect to each of its 
Lagrange multipliers. Inequalities (1.201) are named the conditions of thermodynamic 
stability. 

Defining the specific heat 

1 (dE\ (dS\ (d 2 F\ 

Cv = n [de) V:N = n [de) V:N = ~N {w) VtN (L202) 
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and the isothermal compressibility 

1 fdV\ 1 (dp\ 1 (d 2 G\ 

we see that they are expressed through the second derivatives of thermodynamic po- 
tentials, whence are to be sign-definite. Really, for the specific heat (1.202) we get 

Cv = C ((H 2 ) - (H) 2 ) . (1.204) 



N 

The isothermal compressibility (1.203) can also be found from the expressions 



*-'(*l-'8L-'(k: 



0,N 



Finally, one obtains 

^=^(<^>-<^> 2 ) • (1.205) 
The quantities (1.204) and (1.205) are positively defined, 

C v > 0, K T > . (1.206) 

These are the most often used stability conditions. When at least one of the inequalities 
in (1.206) becomes an equality, this defines the stability boundary. 

1.8 Quasiequilibrium States 

A quasiequilibrium state, or locally equilibrium state, is a statistical state (A) formed 
by the operator averages with a local statistical operator [36,37] (see also [19]). 

To concretize the construction of the latter, introduce the Hamiltonian density 
H(x, t) such that 

Hit) = J H(x,t)dx, (1.207) 
and recall that the number-of-particles operator can be written as 

N(t) = J n(x,t)dx , (1.208) 

where the density-of-particle operator 

n(x, t) = i/j\x, t)tp(x, t) (1.209) 

is defined analogously to (1.80). 

As the local averages (1.130), choose the energy density 

E(x,t) := (H(x,t)} (1.210) 
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and the particle density 
The integration of (1.210), 



p(x,t) := (n(x,t)) . 



(1.211) 



E(t) = J E{x,t)dx , 

gives the internal energy (1.159). The integral of (1.211) yields the average number of 
particles 

N(t) = J p(x,t)dx . 

Note that, according to (1.98), the particle density (1.211) is a diagonal element of the 
reduced density matrix 

pi(x,x,t) = p(x,t) . 

With the given local averages (1.210) and (1.211), the local conditional information 
(1.132) becomes 

I loc (t) = Tr p(t) lnp(t) + X (t) [Tr p(t) - 1] + 
+ J Ai(x, t) [Tr p(t)H(x, t) - E(x, t)] dx+ 
+ / X 2 (x,t)[Tr p(t)n(x,t) — p(x,t)]dx, 



1.212) 



where 



X (t) := In Z loc (t) - 1, Ai(x, t) := f3(x, t), X 2 (x, t) := -f3(x, t)p(x, t) . 

The function fl(x, t) is the local inverse temperature connected with the local absolute 
temperature Q(x,t) by the relation 

P(x,t)Q(x,t) = 1 , 

and p(x, t) is the local chemical potential. 

The principle of minimal information for (1.212) yields the local statistical operator 



Pic 



(t) = 1 exp (- f p(x, t) \h(x, t) - p(x, t)n(x, t)} dx\ , (1.213) 

Zl„ r \t) I J J J 



Hoc{t) 

with the local partition function 



Z loc (t) =Trexp|-y (3(x,t) [H(x,t) - p(x,t)h(x,t) 

Introduce the quasiequilibrium potential Q(t) requiring that 

Z loc {t) = e-^ , 

which gives the definition 

Q(t) :=-\nZ loc (t) = 



dx 



(1.214) 
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= -In Tr exp^- J P(x,t) [6(x,t) - n(x,t)n(x,t)] dx } . (1.215) 
One may also define the local grand potential Q(x,t) by the equation 

Z loc (t) :=exp|-y P(x,t)Q(x,t)dx} . (1.216) 

From (1.215) and (1.216) it follows that 

Q(t) = J p(x,t)Q(x,t)dx . (1.217) 

We need to find how the local grand potential changes when varying the local 
absolute temperature Q(x, t) and the local chemical potential fi(x,t). To this end, we, 
first, obtain from (1.215) the variational derivatives 

5Q ® - E(x,t) - /i{x,t) P (x,t), y\ = -P(x,t)p(x,t). 



5f3(x,t) 5p,(x,t) 
Then, we find the same derivatives from (1.217), 

iQ{t) -IKx,l)+/J(x,l)^M, #1 =/>(*, I)* 1 ** *> 



8f3(x,t) ' ' 5(3(x,t) ' 5fi(x,t) ' 5fi(x,t) 

Comparing the derivatives with respect to P(x,t), we get 

sn(x,t) 



6p(x,t) 
Using the relation 



Q(x, t) [E(x, t) - fi(x, t)p(x, t) - Q(x, t)} 



MM 2 , 5Q(x,t) 

sp(x,t) {x ' } se(x,t) ' 



P(x, t) {tt(x, t) + n{x, t)p(x, t) - E(x, t)] . 



we have 

6Q(x,t) 
8G(x,t) 

Analogously to the equilibrium case (1.179), we may define the local thermodynamic 
entropy 

Substituting to the right-hand side of (1.218) the variational derivative found above, 
we obtain 

S( x , t) = P{x, t) [E(x, t) - fi(x, t)p(x, t) - tt(x, t)] . (1.219) 
Comparing the derivatives of Q(t) with respect to p,(x,t), we find 
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In compliance with (1.218) and (1.220), the variational equation 



6Sl(x,t) 



S(x, t)8Q(x, t) — p(x, t)8p(x, t) 



(1.221) 



holds, named the Gibbs-Duhem relation. 
The local pressure is defined as 



p(x, t) := —Q(x, t) . 



(1.222) 



The variation of (1.219) gives two other useful relations 



P(x,t) 



5S(x,t) 
5E(x, t) ' 



H(x,t) 



e(x,t) 



5S(x,t) 
5p(x,t) 



which show how the local inverse temperature and the local chemical potential are 
connected with the variation of the local thermodynamic entropy. 

In addition to the latter, we can introduce the quasiequilibrium information entropy 



where Q(t) is defined in (1.215). Taking into account (1.217) and (1.219), we obtain 
the relation 



between the quasiequilibrium information entropy (1.223) and the local thermodynamic 
entropy (1.218). 

The local statistical operator (1.213), generally, does not satisfy the Liouville equa- 
tion (1.85), thus, its time evolution cannot be described by a simple unitary trans- 
formation, like (1.86). Therefore, the quasiequilibrium information entropy (1.223), 
contrary to the information entropy (1.119), depends on time. This permits one to use 
the quasiequilibrium ensemble, with the local statistical operator (1.213), as a repre- 
sentation of open, or non-isolated systems whose time evolution cannot be described 
by the Liouville equation. 

When applying the local statistical operator to quasi-isolated systems whose time 
evolution obeys the Liouville equation, one can correct the time behaviour of this op- 
erator by introducing one more Lagrange multiplier, as is described in Section 1.6. 
One may also employ the local statistical operator as a starting approximation in 
constructing perturbation theory for the Liouville equation. Then, the higher-order 
approximations of perturbation theory will automatically correct the starting approx- 
imation adjusting its time behaviour to the Liouville law. Finally, the local statistical 
operator pi oc (0) can be used as an initial condition for the Liouville equation. 



S(t) := -Tr p loc (t) In p loc {t) . 



(1.223) 



With the local statistical operator (1.213), this gives 



S(t) = J (3(x, t) [E(x, t) - p(x, t)p(x, t)] dx - Q(t) , 




(1.224) 
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1.9 Evolution Equations 

The statistical state (A(t)) which is not time-invariant is a nonequilibrium state. The 
statistical operator for which (1.154) is not valid is a nonequilibrium statistical operator. 
The time evolution of the latter is given by the Liouville equation (1.85). 

However, what one finally needs is not a statistical operator itself, but the averages 
of operators associated with observable quantities. For an arbitrary operator A, the 
statistical averages, as defined in (1.87), satisfy the property 

Tr p(t)A(0) = Tr p(0)A(t) , (1.225) 

where A(t) is the Heisenberg representation of an operator A with the time evolution 
(1.81). Therefore, the equation of motion for the average (1.225) can be derived by 
using either the Liouville equation for the statistical operator p(t) or the Heisenberg 
equation for the operator A(t). 

Keeping in mind the general case, when the number of particles in the system is 
not fixed, one defines the generalized Hamiltonian 

H(t) := H(t) - pN(t) , (1.226) 

which can also be called the grand Hamiltonian, since this form is preferable for dealing 
with the grand ensemble. The time evolution is assumed to be governed by (1.226). 

In order to simplify formulas, it is customary to pass to the system of units in which 
the repeatedly appearing Planck constant is set to unity, H — 1. 

In this notation, the Liouville equation for the statistical operator becomes 

l j t p(t) = [H(t),m] , (1-227) 
and the Heisenberg equation for an operator A(t) is 

ij t A(t)=i^A(t) + [A(t),H(t)} , (1.228) 

where the partial time derivative is defined in (1.83). The evolution equations governed 
by a Hamiltonian define the Hamiltonian dynamics. 

The time evolution of a statistical average (A(t)), given by one of the forms in 
(1.225), can be presented by two equations. Using the left-hand side of (1.225) and 
(1.227), one has 

ij t {A{t))=Tr [H(t),p(t)]A(0). (1.229) 
From the right-hand side of (1.225), together with (1.228), one obtains 

l Jt {A[t)) = {t Ft A{t)) + <l A ^ H ^ ■ ( L23 °) 
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These two equations are completely equivalent, and the choice of one of them is a 
matter of taste and convenience [38-42]. According to our own experience, the evolution 
equation (1.230) is usually easier to treat. 

To give examples of evolution equations, we need to concretize the grand Hamil- 
tonian (1.226). For a nonrelativistic system of particles with a two-body interaction 
potential V(x,x'), the generally accepted Hamiltonian is 

H(t) := J VOM) k {V) + U(x, t) - n ifj(x,t)dx+ 

+ I J ^t( X; t)^(x', t)V(x, x')il>(x', t)i>(x, t)dxdx' , (1.231) 

where K(SJ) is a kinetic-energy operator with V:= d/d r, and r:= {r Q | a = 1,2,3} 
being the real-space vector; the multivariable x — { r , . . . , } is assumed to include r ; 
U(x,t) is an external potential; and /i, chemical potential. The field operators satisfy 
the commutation relations 



ip(x,t)i/j*(x',t) =5(x-x f ), 



(x,t)^(x',t)] T = 0, 



1.232) 



which follows from (1.71) and the Heisenberg evolution (1.75), the upper sign corre- 
sponding to Bose-, and lower, to Fermi statistics. Emphasize that (1.232) holds only 
for the field operators at coinciding times. 

From the Heisenberg equation for an annihilation operator 



i^(x,t) = [i;(x,t),H(t)] , 



employing (1.231) and (1.232), one gets 

K(W) +U(x,t) -/I 



■ 9 s 
%—^{x,t) 



+ J §(x,x')^(x',t)^(x',t)^(x,t)dx' , 
where there appears the symmetrized interaction potential 



(1.233) 



1.234) 



having the property 



$(x, x') := - [V(x, x') + V(x', x)} , 



x') = &(x',x) . 



(1.235) 



Note that the same equation (1.234) may be obtained from (1.231) by the variational 
procedure 

i^(x,t) = - 6 , H - . (1.236) 



dt 



6^{x,t) ' 
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The time evolution of the creation operator ip^(x, t) is given by the Hermitian conjugate 
of (1.233)-(1.236). 

To write the equation of motion for the density-of-particles operator (1.209), take 
the kinetic-energy operator as 



K(V) ■= - 



2m 



where m is a particle mass. Introduce the density-of-current operator 

i 



J (x, t) := -— \ ^(x, t) V il>(x, t) - 



V^ f (x,t) 



i/j(x,t) 



(1.237) 



1.238) 



The evolution equation for the density-of-particles operator can be reduced to the 
form 

^n(x,t)+dxv7(x,t) = 0, (1.239) 

called the continuity equation. Here div j=Vj- 

Using again (1.234), together with its Hermitian conjugate, we may derive the 

equations for the components j a (x,t) of the density-of-current operator j (x,t) = 
{j a (x : t)\ a = 1,2,3} given in (1.238). For simplicity, consider the case of no external 
potentials, U (x, t) = 0. Then we find 



d - i d 



dt 



2m 2 p drp 



d 
dr a 



^(x,t) 



d 

drp 



i/>(x,t)+ 





' d 


+ 




pr p 






' d 


m J 




dr. 



-^(x,t) 



d ,. . ld 2 n(x,t) 
2 or a or(i 



dr, 



ij\x, t)ij\x', t)ip(x', t)ip(x, t)dx' . 



(1.240) 



The density-of-energy operator is defined as 

1 



H(x, t) := ^- V il>\x, <)• V il>(x, t)+ 
2m 



1 



4- - / V(x,x')^(x,t)^(x' \t)ip(x' \t)ip(x,t)dx' 



(1.241) 



This is a particular case of the Hamiltonian density (1.207). The evolution equation 
for (1.241) is 



V 2 ^(x, t)- V 4>(x, t)- V ^(x, t) ■ V>(i, t) 



+ 



+ 



i J <f>(x,x')^(x,t)j(x,t)ij(x',t)dx' 
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ij\x', t) j (x, t)i;(x', t) + ^(x, t) j (x', t)il>(x, t) 



dx' . (1.242) 



The evolution equations (1.239), (1.240), and (1.242) are called the local conservation 
laws. Accomplishing statistical averaging of these laws, one derives the equations of 
motion for observables quantities as in (1.230). 

Another way of describing statistical states is based on the fact that these states are 
completely defined by the reduced density matrices with elements (1.96). For example, 
the average number of particles (1.97) involves only 1-matrix. The internal energy 
(1.159), with the Hamiltonian H(t) given by (1.226) and (1.231), that is 

E(t) = {H(t)) = {H{t)+»N) , (1.243) 

can be expressed through 1- and 2-matrices, 

E(t) — J K(x 1 ,x 2 ,t)p 1 (x 2 ,x 1 ,t)dx 1 dx 2 + 

+ \j V ( Xl > x ^P 2 ( X2 > Xl ' X2 ' x i, t )dx 1 dx 2 , (1.244) 
where the notation for the one-particle kernel 

K(xi,x 2 ,t) := 5(xi - x 2 )K 1 (x 2 ,t) 

and the one-particle reduced Hamiltonian 

K^x, t) := K(V) + U(x, t) (1.245) 

are introduced. When the number of particles, N, is fixed, then, using relation (1.101), 
one can cast (1.244) into the form 

E(t) = J 6(xi- yi)5(x 2 - V2) K ^^ jy- P2(xi, x 2 , yi, V2, t)dx x dx 2 dy x dy 2 , (1.246) 

with the two-particle reduced Hamiltonian 

K 2 ( Xl ,x 2 ,t) := K^xut) + K^t) + (N - l)V(x u x 2 ) . (1.247) 

To find the internal energy (1.246), one needs to know only 2-matrix. 

Generally, all averages corresponding to observable quantities can be expressed 
through several, usually a few, low-order density matrices. Therefore, one has to find 
out the evolution equations for the latter. These equations can be derived by using 
(1.234). Thus, for 1-matrix we have 

d 

i Tr:Pi(xi,x 2 ,t) = [Ki(xi,t) - K 1 (x 2 , t)\ pi{x u x 2 , t)+ 
at 

+ j [$(xi,x 3 ) -®(x 2 ,x 3 )]p 2 (x 3 , xi,x 3 , x 2 ,t)dx 3 . (1.248) 
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For 2-matrix, we obtain 



d 

i-Qj.P2(xi, x 2 , yi, V2, t) = [H(x u x 2 , t) - H(y u y 2 , t)] p 2 {x u x 2 , yi, y 2 , t) + 

+ J [$(xi,x 2 , x 3 ) - $(yi, y 2 , x 3 )] p 3 (x 3 , x u x 2 , x 3 , y u y 2 , t)dx 3 , (1.249) 
where the two-particle Hamiltonian 

H(x u x 2 , t) := K^xu t) + K^x-2, t) + $(x u x 2 ) , (1.250) 
and the three-particle interaction 

$(xi, x 2 , x 3 ) := x 3 ) + $(x 2 , x 3 ) 

are introduced. 

If the kinetic-energy operator has the form of (1.237), then, using the property 

v? - v\ = (vi + v 2 ) (vi - v 2 ) , 

the 1-matrix equation (1.248) can be written as 

i^.Pi(xi,x 2 ,t) = -i (Vi + V2) j (x 1 ,x 2 ,t) + [U(xi,t) - U(x 2 ,t)]p 1 (x 1 ,x 2 ,t) + 



+ J [^(x 1 ,x 3 ) -$(x 2 ,x 3 )]p 2 (x 3 , Xi,x 3 , x 2 ,t)dx 3 , 
where the density-of-current matrix is defined by 



(1.251) 



j (xi,x 2 ,t) 



Vi - V; 

2mi 



pi(x 1 ,x 2 ,t) . 



(1.252) 



The relation between (1.252) and the density-of-current operator (1.238) reads 

(j(x,t)) = lim j (x,x',t) . 

x'— 

The evolution equation for the matrix in (1.252) is 
d - 

1— j (x 1 ,x 2 ,t) = [Kifaut) - K 1 (x 2 ,t)\ j (x 1 ,x 2 ,t)- 



1 

2m 



Vi U{x u t)+ V2 U(x 2 ,t) 



p 1 (x 1 ,x 2 ,t)- 



1 

2m 



Vi ^(x!,x 3 )+ V2 ®(x 2 ,x 3 ) 



p 2 (x 3 ,x 1 ,x 3 ,x 2 ,t)dx 3 - 



49 



J [<&(xi, x 3 ) - $(x 2 , x 3 )] (Vi - V2) P2(^3, X!, x 3 , x 2 , t)dx 3 , (1.253) 



% 

2m 

where the 2-matrix satisfies (1.249). 

The evolution equations for matrices are to be supplemented with initial and bound- 
ary conditions. However, the system of these equations is not closed, since the equation 

for pi or j involves p 2 , the equation for p 2 involves p 3 . Generally, an equation for p n 
always contains p„+i, of course, if the interaction potential is not zero. In this way, re- 
duced density matrices satisfy an infinite system of integrodifferential equations, which 
is called the hierarchical chain of equations. To obtain a closed set of equations one 
has to break this chain by approximately expressing one of p n through the lower-order 
density matrices. Then one gets a closed system of n — 1 equations for p 1 , p 2 , . . . , p n -i- 
This kind of approximation is called the decoupling procedure. It is important that any 
approximation used would preserve the general properties of density matrices following 
from their definition (1.96). This is the conjugation property 

p* n (x n ,y n ,t)=p n (y n ,x n ,t); (1.254) 

the symmetry property 

P(x i ,x j )p n (x n ,y n ,t) = ±p n (y n ,x n ,t) , (1.255) 

where ; the non-negativeness of the diagonal elements 

Pn(x n ,x n ) >0; (1.256) 

and the normalization condition (1.103). For the density-of-current matrix, given by 
(1.252), in compliance with (1.254), one has 

— >* — > 

j (x l ,x 2 ,t) =j (x^x^t) . 

Note that in some cases the symmetry property (1.255) can be broken giving, nev- 
ertheless, a good approximation [43-46]. For example, the widely used Hartree approx- 
imation 

P2(xi,x 2 ,yi,y2,t) = pi(x 1 ,y 1 ,t)p 1 (x 2 ,y 2 ,t) , (1.257) 

obviously does not obey (1.255). The property (1.255) holds for the Hartree-Fock 
approximation 

p 2 (x u x 2 ,y u y 2 ,t) = pi(x 1 ,y 1 ,t)p 1 (x 2 ,y 2 ,t) ± p 1 (xi,y 2 ,t)p 1 (x 2 ,yi,t) . (1.258) 

Substituting a decoupling, like (1.257) or (1.258), into (1.248), one obtains a closed 
equation for 1-matrix, 

d 

i-^Pi(x 1 ,x 2 ,t) = [K 1 (x 1 ,t) - K 1 (x 2 ,t)]p 1 (x 1 ,x 2 ,t) + 
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+ J [<&(zi,z 3 ) - ®(x 2 ,x 3 )] [pi(x 1: x 2: t)p 1 (x 3 ,x 3 ,t)+ 

+ep 1 (x 1 ,x 3 ,t)p 1 (x 3 ,x 2 ,t)]dx 3 , (1.259) 

in which e = corresponds to the Hartree approximation and e — ±1, to the Hartree- 
Fock approximation. Initial and boundary conditions are supposed to be given for 
equation (1.259) whose solution has also to satisfy the conditions 

p* 1 (x 1 ,x 2 ,t) = pi(x 2 ,X!,t) , pi(x,x,t) > , 

^ y p 1 (x,x,t)dx = p , (1.260) 

where p is an average density of particles. 

Actually, equation (1.259) is a set of two equations, since 1-matrix is a complex 
function 

Pi(xi, x 2 , t) := a(xi, x 2 , t) + ib(xi, x 2 , t); (1.261) 

here a and b are, respectively, the real and imaginary parts of p 1 . Owing to (1.260), 
we have 

a(x 2 , xi, t) = a(x 1 ,x 2 , t), b(x 2 , Xi, t) = —b(x 1 ,x 2 , t), 
a(x,x,t) > 0, b(x,x,t) =0, —J a(x,x,t)dx = p . 

Then, from (1.259) and (1.261) we get 
d 

—a(x 1 ,x 2 ,t) = [K 1 (x 1 ,t) - K 1 (x 2 ,t)]b(x 1 ,x 2 ,t)+ 

+ J [$(xi,x 3 ) - &(x2,x 3 )] {a(x 3 ,x 3 ,t)b(xi,x 2 ,t)+ 

+e [a(x 1 ,x 3 ,t)b(x 3 ,x 2 ,t) + a(x 3 , x 2 ,t)b(x!, x 3 ,t)]} dx 3 , 
d 

—b(x 1 ,x 2 ,t) = [Ki(x 2 ,t) - K 1 (x 1 ,t)]a(x 1 ,x 2 ,t)+ 

+ J [&(x 2 ,x 3 ) - ®(x 1: x 3 )] {a(x 1 ,x 2 ,t)a(x 3 ,x 3 ,t)+ 

+e [a(x 1 ,x 3 ,t)a(x 3 ,x 2 ,t) + b(x± : x 3: t)b(x 3: x 2 , t)]} dx 3 . 

These equations describe the dynamics of all, diagonal as well as off-diagonal, elements 
of 1-matrix. 

In addition to the Hartree approximation (1.257) and Hartree-Fock approximation 
(1.258), another widely used decoupling is the Hartree-Fock-Bogolubov approximation 
[39] 

p 2 (x 1 ,x 2 ,y 1 ,y 2 ,t) = p 1 (x 1 ,y 1 ,t)p 1 (x 2 ,y 2 ,t)± 

± Pi (xi, 2/2, t) Pl (x 2 , Vl ,t) + (ij(x 1 )^(x 2 ))(^(y 2 )^(y 1 )) (1.262) 

containing the anomalous averages (ipip) and {^rtjr). In this case, to find a closed set of 
evolution equations for 1-matrix, one needs to write down the evolution equations for 
the anomalous averages. Group-theoretical analysis of the Hartree-Fock-Bogolubov 
approximation reveals a rich variety of admissible physical states [47]. 
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1.10 Coherent States 

The notion of coherent states appears in many applications, therefore it is worth pre- 
senting here the general mathematical definition and the main properties of these states 
[48-51]. 

Coherent states are microstates that are normalized eigenvectors of the annihilation 
operator. Taking the latter in the form (1.54), we have, by definition, 

ip(x)h v := i](x)h v , (1.263) 

where r)(x) is a coherent field, and the coherent state 

h v = [ Vn (x n )] e H (1.264) 

pertains to the Hilbert space H given by (1.51). 

The action of the annihilation operator (1.54) on the vector (1.264) yields 



i)(x)h n = \Vn+ lr] n+1 (x n ,x) 
Substituting this into (1.263), we have 

Vn + lrj n+1 (x n ,x) = i](x)r] n (x n ) . 
The solution to this recursion relation reads 

Vn (x n ) = ^= n v(x k ) , (1.265) 

where r] G C and rj(x) is an arbitrary complex function of a Hilbert space with the 
scalar product 

The scalar product of coherent states is 

n=0 n=0 n - 

which can be reduced to 

From the normalization condition 



h+h„ = 1 



we find 



|?7o| = exp|-^(?7,77)| . (1.266) 
Thus, coherent states are defined by (1.264), (1.265), and (1.266). 
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For the scalar product of coherent states, within a constant phase factor, we have 
K h 6 = e*P {-|(77, 77) + , (1-267) 

which shows that different states are not orthogonal to each other, since 

< \h+hz\ < 1 . 

Introducing the corresponding functional integral with a differential measure Vrj, one 
can get the resolution of unity 



h v h+V V = 1 



(1.268) 



The formal expression (1.268) is to be interpreted in the sense of weak equality, that 
is, as equality of arbitrary matrix elements of the indicated expression. 

Coherent states form a basis but it is overcomplete because of the absence of or- 
thogonality. Also, coherent states are eigenvectors of only the annihilation operator 
(1.54), but not of the creation operator (1.55), whose action on (1.264) gives 



fthr, = [Vn7] n _ l (x n 1 )5(x n -x) 
However, there is a relation 

h+ft(x) = h+r]*(x) 

following from the Hermitian conjugation of (1.263). 
Considering the double action 

^(x)^(x')h v = 77(2) 77(2/)^? , 
we see that, on the class of coherent states, annihilation operators commute, 



1.269) 



tjj(x),ijj(x') 



h v = . 



But there are no simple commutators involving creation operators, which is clear from 
the actions 



^(x)^(x')h v 



n(n - l)i] n _ 2 {x n )5(x n -i - x')8(x n - x) 



■4)(x)-4)\x')h ri = [(n + l)r) n (x n )5(x - x')] , 
^ f (x')^(a;)/ir, = nr] n (x n ~ 1 ,x)5(x n - x') . 
Notice that function (1.265) is symmetric with respect to permutation of any two 

variables Ob ^ j Ob j Ob j 

P{x u x j )r jn (x n ) = r] n (x n ) . 
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Consequently, rj n (x n ) G and h v G 7i + , where 7i+ and 7i+ are defined in Section 
1.3. In the same way as in (1.69), we may introduce the Bose field operators ip + and 
ifj\_, which have, on the class of coherent states, simple commutation relations given in 
(1.71). 

To introduce coherent states for the Fermi field operators, according to (1.69), 
requires that the function (1.265) would be antisymmetric. This is achievable if the 
coherent field r](x) is not a complex function, as in the case of Bose statistics, but an 
element of the Grassmann algebra, satisfying the relation 

r) F (x)r) F (x') + r) F (x')r) F (x) = , 

in which the index F refers to the Fermi case. If r](x) is a Grassmann variable, then the 
function (1.265) is antisymmetric with respect to permutation of variables 
Hence, rj n (x n ) G H~, and h v G H-. Therefore, one can define the Fermi field operators 
if)- and as in (1.69), which have, on the class of coherent states, the standard 
anticommutation relations, as in (1.71). 

In this way, coherent states can be introduces as the normalized eigenvalues of Bose, 
as well as Fermi operators, 

ip(x)h v := rj(x)h v , (1.270) 

where ip(x) is either ip + or tp-(x). The difference is that for Bose statistics rj(x) is 
a complex function, while for Fermi statistics r)(x) is an element of the Grassmann 
algebra. These two cases may be combined in one condition for the coherent fields 

[ V (x), V (x')] T = 0, (1.271) 

in which the upper sign, as usual, is related to Bose statistics, and the lower, to Fermi 
statistics. 

Coherent states (1.264) can be obtained from vacuum by means of field operators. 
By analogy with (1.57) and (1.63), we may define n-particle coherent states 

KnVn(x n )} := 4=f / Vn(x n ) f[ ^(x k )\0)dx n , (1.272) 

where r](x) entering into (1.265), depending on the type of statistics considered, is 
either a complex function or a Grassmann variable, as is discussed above, and i]y{x) 
is, respectively, either a Bose creation operator, or Fermi creation operator. With the 
form (1.265), we can present (1.272) as 

[6mnVn(x n )} = ^ (/ ^ {x^dxj |0> . 

Then, for the coherent state (1.264) we get 

h v = jr^([^(x)v(x)dx) n \o). 
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As is easy to notice, the latter expression is a series expansion of an exponential 

h v = ^o ex p|y ip^ (x)r)(x)dx^ |0) . 

Taking onto account the normalization (1.266), we come to 

/i„ = exp{-l(j7,J7) + J ^(xHx)dx}\0) . (1.273) 

This is a general formula for constructing coherent states from vacuum. 

To find the time evolution of coherent states, let us generalize the definition (1.263) 
to the form 

ip(x, t)h v = rj(x, t)h v , (1.274) 
in which i/)(x,t) is given by (1.75). Whence, the time-dependent coherent state is 

h v (t) := U(t)h v , (1.275) 

for which 

^(x,0)h v (t) = r](x,t)h v (t) . 
The Hermitian conjugate to (1.274) is 

h^(x,t) = tf V *(x,t) . 

The evolution equation for (1.275) follows from the equation 

ijU(t) = HU{t) 
for the unitary operator, which yields 

ijh n (t) = Hh n (t) . (1.276) 

The scalar product of time-dependent coherent states 

h+(t)ht(t) = h+ht = rj^ e^ 

does not depend on time. 

Coherent states permit us to write simple expressions for the matrix elements of the 
operators from the algebra of observables. Thus, for the number-of-particle operator 
(1.208), using 

h^(x,t^(x',t)h^ = ^V^VOM)^',*), 

we get 

h+N(t)h^ = rifoeW f V *(x,t)Z(x,t)dx . 
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The diagonal elements of this is 

h+N(t)h v = J \r]*(x,t)\ 2 dx , (1.277) 

where the normalization (1.266) is invoked. For the Hamiltonian (1.226), the diagonal 
element reads 



h+H{t)h n = J rf{x,t) 



K(v)+U(x,t) - ii 



r)(x, t)dx+ 



+ \rj(x,t)\ 2 V(x,x')\r](x',t)\ 2 dxdx' . (1.278) 

The evolution equation for a coherent field r)(x, t) can be derived from equation 
(1.234) for the field operator and the use of the equality 

h+i^(x, t)h v = i^ V (x, t) . (1.279) 

In this way we obtain the coherent-field equation 

^K(V) + U(x,t)-/2 + J<f>(x,x')\ri(xr,t)\ 2 dx'}ri(x,t) = i^q{x,t) . (1.280) 

The latter has the structure of a Schrodinger equation with an effective coherent po- 
tential 

U eff (x,t) := U(x, t) - ii + / $(x,x')\r)(x',t)\ 2 dx' . 

Since the part of (1.280), containing interparticle interactions is nonlinear with respect 
to the coherent field r)(x,t), the equation (1.280) is a nonlinear Schrodinger equation. 
This equation defines the coherent field up to a phase factor. 

It is important to notice that nontrivial coherent states exist not for any interaction 
potential, but only for those that are integrable, so that 



<&(x,x')dx' 



< oo 



If it is not so, then the corresponding term in (1.280) diverges and the coherent-field 
equation has only a trivial solution rj(x,t) = 0. 
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2 Green's Functions 



The method of Green's functions is one of the most powerful techniques in quantum 
statistical physics [52-57] and quantum chemistry [58]. 

If L is a linear operator, and f(x) is a given function, then a solution of Lu = f 
can be expressed in the form u(x) = J G(x, £)f(£)d£, where G satisfies the equation 
LG(x,C,) = 5(x — £). The function (2(1,2) which serves as the kernel of the above 
integral operator, is called a Green's Function. This, because it was George Green 
(1793 — 1841) who introduced these functions into the study of linear equations. Certain 
Green's functions are also called propagators because as follows from the above integral 
we can think of G as transfering an L-modulated effect of a "charge" /(£) at £ to the 
point x. 

In recent years, the use of Green's functions has become widespread in theoretical 
physics and especially in the study of the iV-body problem in quantum mechanics. As 
we shall see below they occur in many different forms. A basic mode of classifying them 
is as one-particle, iu>o-particle, . . ., and n-particle propagators. A reduced density 
operator can be defined as an initial value of a corresponding type of a Green's function 
and so the latter are also restricted by the necessity of being N-representable. This 
fact was often neglected by early Green's function enthusiasts. For a statistical system 
a hierarchy of equations, involving n— and (n+1)— particle Green's function from n = 
to n = N, must be satisfied. Only thus can iV-representability of the Green's functions 
be guaranteed. For large N this is manifestly impossible, so in practice the hierarchy 
of equations is usually truncated by an ansatz at the level of the first, second or third 
order Green's functions. But this always leaves some doubt as to the extent of possible 
errors introduced by the ansatz. 

An even more difficult obstacle in discussing the iV-body problem is that it cannot 
be solved exactly if AT > 1! So to make any progress, some form of approximation 
must be adopted. For 2 or 3 particles, the approach which Hylleraas used to explain 
the spectrum of atomic helium can sometimes give very accurate results. However, 
for larger N one normally starts from a tractable initial approximation attempting 
to improve this by successive perturbations. During recent decades physicists have 
developed increasingly sophisticated methods for obtaining approximations to Green's 
functions, and for the application of perturbation techniques to them. 

It has been our aim in this chapter to provide the reader with a comprehensive 
overview of the extensive resulting literature about Green's functions. We have not 
followed any one previous treatment. Rather, we have sought to gather together what 
we regard as the most important and most generally applicable aspects of the Green's 
function approach to the quantum mechanical iV-body problem. It may be that, for 
some readers, our condensed treament will bring to mind the well-known canard: "in 
order to do theoretical physics all you need is a knowledge of perturbation theory and 
a mastery of the Greek alphabet" and they will feel submerged under formalism! If so, 
we apologize, and only plead that we were unable, otherwise, to present in reasonable 
length what we consider to be the state-of-the-art for Green's functions. 
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2.1 Principal Definitions. 



We shall adopt the widespread convention of using Planck's constant divided by 2n 
and Boltzmann's constant as units; so h/2ir = 1, ks — 1. 

In order to have a compact notation we shall often denote the particle coordinates 
(xi,ti) by the index i, so that a function of n particles may be written as 

f(\2 ...n):= f(x u ti, x 2 , t 2 ;...; x n , t n ), 

and a differential for n particles becomes 

d{\2 . . .n) := dx\dt\dx 2 dt 2 . . .dx n dt n . 

Thus, for the field operator we shall often use 

:=^(a:i,0- 

The distinction between Bose and Fermi statistics is basic to quantum mechanics and is 
encoded in the field operator. We have adopted the convention that when the symbols 
± or =F appear in our formulas the upper sign refers to bosons. 

It is a basic metaphysical assumption of most of us that cause precedes effect in 
time. In consequence a so-called chronological or time-ordering operator, T, acting on 
field operators plays an important role in the theory of Green's functions. 

TM1UU2)-1 ^W^ 2 )' tl2>0 (2 1) 

iVW w -j±^t (2Wl) , tl2<0 ' v- L ) 

where 

tu '■= ti — t 2 = —t 2 \. (2-2) 

Similarly, when T acts on the products ip(l)ip(2) and ^(1)^(2), it arranges the 
order of the field operators so that time increases from right to left. If a permutation 
of the operators is neccessary, the product is multiplied by +1 or —1 according as Bose 
or Fermi statistics is relevant. 

Invoking the unit-step function: 

<*>=-{£:>!• 

the time-ordered product (2.1) can be written as 

f ^(1)^(2) = (t 12 ) V (1) V f (2) ± O (-* 12 ) ^ (2) i> (1) . 

The time-ordered products of any number of field operators can be defined in a similar 
manner. 

We now define five different types of Green's functions commonly used in many- 
body theory: 
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• causal Green's function or propagator 

G(12) := -i < Ty>(l)V f (2) >; 



(2.3) 



• retarded Green's function 



G ret (12) :=-iG(t 12 ) < V(l)^ j (2) >; 



(2.4) 



• advanced Green's function 



G adv {12) :=-20(t 21 )<^(2)^(l)>; 



(2.5) 



• commutator and anticommutator retarded Green's function 



:= -O (f 12 ) < [^(1),^(2)] >; 



(2.6) 



• commutator and anticommutator advanced Green's function 



G^(12) :=-i0(t 21 )< ^(l),^(2)j T >. 



(2.7) 



Note that the average < • > may , depending on the context, be with respect to a 
pure state or with respect to an ensemble which may be temperature dependent and 
which may be in equilibrium or not. 

The properties of these different Green's functions are similar to one another and 
there are fairly obvious relations among them so that in most situations it is as conve- 
nient to work with one of them as another. However, in our opinion the causal Green's 
function has a slight edge for the following reasons: 

1. As a linear combination of the retarded and advanced Green's functions, the 
causal Green's function carries more information than either of those separately: 



2. Dealing with the causal Green's function is rather simpler than with the commu- 
tator or anticommutator Green's functions; 

3. The causal Green's function, as we have defined it is accepted and used, almost 
universally, in both quantum statistical mechanics and quantum field theory. 

Notice, further, that in the usual set-up in Hilbert space, operators act on the left 
of the symbols denoting states of the system. Thus when the operator is a product we 
interpret its action to mean that the factor on the right acts first. The chronological 
operator, T, was defined with this in mind ensuring that the factor which acts first on 



G (12) = G ret (12) ± Gad v (12) ; 



(2.8) 
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the wave vector is at the earlier time. This corresponds to our intuitive feeling about 
causality, justifying the term "causal" applied to (2.3) 

It turns out that G (12) is closely associated with the time-dependent 1-matrix. In 
fact, 

Pi(xi,a; 2 ,t) = ±iliinG(12), (2.9) 

where the limit is taken subject to t 2 > t\. 

Recall that here and elsewhere when ± or =f occur indicating Fermi or Bose stastis- 
tics, the upper sign refers to Bose statistics. 

The diagonal elements of the 1-matrix give the value of particle density and are 
therefore positive, that is 

pi (x,x,t) = p (x,t) . 
From this we obtain the following condition on G(12), 

iilim lim G (12) > 0, (2.10) 

where the limit is taken subject to the restriction t 2 > t\. 

All of the Green's functions defined above are one-particle functions. In an anal- 
ogous manner we may define a two-particle causal Green's function, or two-particle 
propagator. 

G 2 (1234) :=-<fij} (1) tp (2) ift (3) ift (4) > . (2.11) 

The definition (2.1) for the chronological operator is generalized to an arbitrary 
finite product of field operators as follows: the factors are permuted so that time 
increases from right to left and the appropriate sign is determined by decomposing the 
required permutation into a succession of transpositions of neighbouring factors. It 
then follows that (2.11) has the following symmetry properties, 

G 2 (1234) = ±G 2 (2134) = ±G 2 (1243) . 

We again find that there is a simple relation between the 2-particle propagator and 
the time-dependent 2-matrix, namely 

p 2 (x 1 ,x 2 ,x 3 ,x 4 ,t) = T lim G 2 (1234), (2.12) 

where the limit is taken subject to the restriction t 4 > t 3 > t\ > t 2 . 
Similarly, we define the n-particle propagator as 

G n (12 . . . 2n) := (-i) n < Tty (1) . . . -0 (n) ^ (n + 1) . . . ^ (2n) > . (2.13) 

Again, we can recover the n-matrix from (2.13). It turns out that 

Pn (X 1: . . . , X 2n , t) = (±l)(3n-l)n/2 -n y ^ (1 . . . 2 n) , (2.14) 

where the limit is subject to the restriction 

^2n > t 2n -i > . . . > t n+ i > t\ > t 2 > . . . > t n . 
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In all of the Green's functions defined so far there are the same number of creation 
and annihilation operators. This is appropriate for situations in which the number 
of particles is conserved for in such cases the expected value for operators involving 
unequal numbers of creation an annihilation operators will be zero. However, if the 
number of particles is not conserved then such averages may be different from zero. 
Indeed, this can happen for systems with a Bose condensate and for superconductors. 
In these cases averages of operators such as ijj(l)ip(2), ^(l)^ (2), ^ {l)ip(2)ip{?>) may 
differ from zero. Such averages, which can differ from zero only if the number of 
particles is not conserved are often called anomalous averages and the Green's functions 
defined by means of such averages are called anomalous Green's functions. The most 
commonly used are certainly 

F(12) := -% < TV(1)V(2) >, 

F + (12) := i < 1^(2)^(1) > . (2.15) 

In analogy with the definition of Green's functions by means of the field operators, 
it is possible to employ any other operators to define Green's functions as statistical 
averages of products of operators. All observable quantities can be expressed by means 
of integrals involving Green's functions. For example it follows from (2.9) that the 
density of particles is 

p(s,f)=±KmtonG(12). (2.16) 
Here we have employed an abbreviated notation for the limit to mean 

hm^x 2 -^x 1 ; U->t, (t 2 >h). (2.17) 

(21) 

Further, the density of particle current can be expressed as 

JOM) = g?f( v i - V 2 )G(12) (2.18) 

where the limit is interpreted by means of (2.17). 

Any 1 -particle operator, Ai(t), can be expressed by means of the field operators 
and an operator- valued function A 1 (x, t) as 

Ai(t) = J ^(x,t)A 1 (x,t)ip(x,t)dx, 

and its average or expected value can be written as 

< A 1 >= ±i J \imA 1 (x 1 ,t)G(12)dx 1 . (2.19) 

Similarly, the average, or expected value, of a 2-particle operator 

A 2 (t) = 7} J ^(x 1 ,t)4)\x 2 ,t)A 2 (x 1 ,X2,t)4j(x2,t)4)(x 1 ,t)dx 1 dx2 
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IS 

< A 2 > 



~ I lim A 2 (x 1 ,x 2 ,t)G 2 (123A)dx 1 dx 2 , (2.20) 

2 J (4321) 



(4321) 

where in this limit: 

x 4 — >xi, x 3 -> x 2 ; U—>t, (U > t 3 > h > ti). (2-21) 

In order to abbreviate notation, we shall use x n , with n as a superscript, to denote 
the n variables (x±,x 2 , . . . ,x n ) and thus dx n := dx\dx 2 ...dx n . Then an n-particle 
operator can be written as 

Mt) = — [ i/>\ Xl ,t) . . . 4>\x n , t)A n {x n , t)^j{x n , t)... 4>( Xl ,t)dx n . 
The average of this operator is given by the formula 

< A n >= ^—J- J Urn A n (x n , t)G n (l . . . 2n)dx n (2.22) 

which involves the n-particle propagator (2.13) and the generalization of the limiting 
procedures (2.17) and (2.21) defined by the rule 

lim := J O^^-. 1 (2.23) 

(2n...l) { tj -> t, t j+1 > tj, Vj. 

Applying this to the number-of-particles operator, N, we find that the average 
number of particles is 

N =< N >= ±i J limG , (12)dxi. (2.24) 
Averaging the Hamiltonian, H, we find 

<H >=±ij lim [K(Vi) + C/(xi, i) - u] G(12)ttei- 

--/ lim Wxi,a:2)G : 2(1234)da;ida;2- (2.25) 

2 7 (4321) 

Using this we obtain the internal energy as 

E=<H > +fiN. (2.26) 

These examples should convince the reader that an exact knowledge of the Green's 
functions for an iV-body system contains complete information about its state whether 
or not it is in a pure or statistical state. However, we also know that for iV > 1 it is 
too much to hope for exactness. So our problem is to develop methods which can give 
us reasonable approximations. 
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2.2 Spectral Representations 



The Green's functions defined in the previous section are appropriate whether or not 
the system being studied is in equilibrium. When it is in equilibrium, as suggested by 
Lehmann, taking the Fourier transform with respect to the time variable leads to a 
form for the Green's functions which is referred to as the spectral representation and 
with which it is much easier to work. The properties of this representation give rise 
to certain advantages of the Green's function method. In particular, singularities in 
the spectral representation of the first order Green's function occur at energy levels of 
quasi-particles. 

For the equilibrium states we can introduce the equilibrium correlation functions: 



C + (x u x 2 ,t) :=< ^(1)^(2) >, 
C4x u x 2 ,t) :=< V f (2)V(l)>, 



in which 



t = t 



12 



to. 



[2.27) 



(2.26 



The right hand side of the equations in (2.27) apparently depends on t\ and t 2 but, 
in fact, only on the difference t. For recall that when the Hamiltonian, H , does not 
depend explicitly on t, i.e. dH/dt = 0, the evolution operator is given by 

U(t) = e~ mt ; p(t),p]=0, 

where p is the iV-particle von Neumann density matrix of the system. By definition, 

< ^(1)^(2) >= Tr [^(1)^(2)] =< iP(x 1 ,t)^(x 2 ,0) > 

since ip(x,t) = U + (t)tp(x, 0)U(t) 

The functions (2.27) enjoy the following properties with respect to complex conju- 
gation. 

C* + (x 1 ,x 2 ,t) = C + (x 2 ,x 1 ,t) Cl(x 1 ,x 2 ,t) = C-(x 2 ,x 1 ,t). (2.29) 

Recall that f3 denotes the inverse temperature. An important relation between the 
correlation functions (2.27) is the equilibrium condition: 

C-(x!,x 2 ,t) = C + (x u x 2 ,t-il3), (2.30) 

which is often referred to as the Kubo-Martin-Schwinger conditon. Its proof follows: 



C_( Xl ,x 2 ,t) =Tr\pij\x 2 ,0)ij( Xl ,t) 



Tr 



ip(x 1 ,t)pip\x 2 ,0) 



Tr 



pe^ H ^(x l7 t)e- pH ^(x 2 , 0)1 = Tr[p^(x l7 t - i(3)^(x 2 , 0)] = C + (x u x 2 , t - i(3) 



since, at thermal equilibrium, p = (exp(— PH))/Tr(exp(— f3H)). 
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The Fourier transforms, K±, of functions (2.27) are given by 

/+OO . Ajjj 

K + (x 1 ,x 2 ,uj)e~ tu ' t — , 
-oo Z7T 

/+oo 
C + (x 1 ,x 2 ,t)e lu ' t dt , , 
-oo 

/+oo A^j 
K4x 1 ,x 2 ,cu)e-™ t — , 
-oo Z7T 

f+OO . 

K_(xi,x 2 ,u)= C^(x 1 ,x 2 ,t)e tuJt dt . 

J — oo 

Then, taking the Fourier transform of (2.30) we obtain 

K_ (x 1 ,x 2 , u) = e~^K + (x ± ,x 2 , u) . (2.31) 
The correlation functions (2.27) at t — satisfy the equation 

C+ ( Xl ,x 2 , 0) =F C_ (x u x 2 , 0) = 5 ( Xl - x 2 ) , (2.32) 

which follows from the commutation relations for the field operators. The Fourier 
transform of (2.32) reads 

/+oo , . Aijj 

^ K + (x 1 ,x 2 ,u;)(lTe^) — = 5(x 1 -x 2 ), (2.33) 

when we take (2.31) into account. 

By means of the correlation functions (2.27) it is possible to rewrite the Green's 
functions of (2.4) and (2.5) as 

G ret (12) = -iG(t)C + (x 1 ,x 2 ,t), G adv (12) = -iG(-t)C-( Xl , x 2 , t), (2.34) 

where the notation (2.28) is used. The Fourier transforms, G + and G_, of these Green's 
functions are defined by the relations: 

G ret (12) = — / G + {x u x 2 ,u)e-^du, 

Z7T J -oo 
\ r+oo 

G adv (12) = — G-(x u x 2 ,tu)e- wt du. (2.35) 

Z7T J-00 

The unit step function, 0, can be expressed as a Fourier integral as follows: 

/+oo p -iut J 
— , (2.36) 
-oo cj±i0 2vr' V ; 

where the convention 

/+oo r+oo 
f(u±i0)du : = lim / / (u ± ie) du (2.37) 
-OO £—++0 J ^QQ 
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is used. 

Recall the well-known formulas for the 5-function: 

/ + OO r\i.) r+oo 

e-"*^- / e wt dt = 2h5(uj). 
-oo 2n J-oo 

To obtain expressions for G + and G_, we invert equations (2.35) making use of (2.34) 
and the convolution theorem for Fourier integrals. This leads us to the spectral repre- 
sentation of the retarded and of the advanced Green's functions 

G+(xi,x 2 , u) = -— / — — alu , 

2n J-oo uj' — uj — lO 

G-(x 1 ,x 2 ,u) = — — ——duo . 2.38 

2n J-oo uj' — uj + lO 

The spectral representation of the causal Green's function is an immediate conse- 
quence of (2.38) and the definition of the Fourier transform of (2.3) as 

/+oo r+oo dm 

G(l2)e tuJt dt, G(12)= G(x 1 ,x 2 ,uj)e- tuJt — . (2.39) 
-oo J —oo Z7T 

According to (2.8) we have 

G (x 1 ,x 2 , uj) = G + (x 1 ,x 2 , uj) ± G_ (xix 2 , uj) ■ (2.40) 

Combining (2.38) and (2.40), we obtain the spectral representation 

r+oo ( i e~ i/3bJ ' \ duj' 

G(x u x 2 , uj) = - K + (x u x 2 , uj') =F — - — , (2.41) 

J-oo \uj' — uj — lO uj ' — uj + lO J 2n 

of the causal Green's function. 

For systems in equilibrium, the averages of operators are expressed by means of 
integrals involving G(xi,x 2 ,u>) using the formula (2.41). For example, the density of 
particles (2.16) is 

/+OO Am 
e^G(x,x,uj) — . (2.42) 
-oo ZTT 

The Green's function at coincident times is defined as 

G(12) | t=0 := t hmG(12). (2.43) 

Substituting (2.41) into (2.42), we obtain 

" +0 ° „ , s _ a .,dcj 



p(x,t) = / K + (x,x,uj)e-^ — . 

J-oo ZTC 



Instead of dealing with the two Fourier transforms K + and K_ which are related 
by (2.31), we introduce the spectral function 

l±n(u) n{uj) 
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Here, 

: = ^Vt (2 - 45) 

is called the i?ose function or the Fermi function according as one takes the minus or 
plus sign. Recall that this is consistent with our general convention that the upper of 
the two signs refers to Bose statistics, the lower to Fermi statistics. 

Using the equality l±n(uj) = n(uj)e l3u) and the spectral function (2.44) the spectral 
representation of the causal Green's function, or propagator, can be expressed in the 
form 

/+oo 
J(xi,x 2 ,ou') 
-oo 

Expressions involving poles can be reorganized by invoking the symbolic identities, 
introduced by Dirac, in which the symbol P indicates that the principal value is to be 
taken. ^ ^ 

:= P— =F iirS(u)), 

tu±iO uj y h 

f(w±iO) f(uj) 
Evaluating (2.46) by means of (2.47), we find that 

G{x u x 2 , u)=P J ^ x ^f) ^L _ » [1 ± 2n{u)\ J( Xl ,x 2 , u). (2.48) 

J-00 UJ — LO ZTT Z 

If J has no poles on the real axis then again using (2.47) we can show that (2.48) 
implies that 

J (xi, x 2 , u) — i[G (x\x 2 , uj + iO) — G (xi, x 2 , uj — iO)} . (2.49) 



l±n(uj') 



UJ' 



UJ 



iO 



n{uj') 



UJ' 



uj + iO 



) ■ 



(2.47) 



2.3 Dispersion Relations 

As follows from their spectral representations (2.38) the Fourier transforms of the 
retarded and advanced Green's functions, G + and G-, have the following asymptotic 
behaviour, when uj — ► ±00 

G+ (x lt x 2 , uj) = -C + (x!,x 2 , 0), 

UJ 

and 

G- (xi, x 2 , uj) = C- (xi,x 2 , 0) . 

UJ 



Thus 



lim G + (xi, x 2 , uj) — lim G-(x±, x 2 , uj) — 0. 
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It follows from (2.34) that the retarded Green's function differs from zero only if t > 
0, and the advanced Green's function only if t < 0. From this and the expression (2.35) 
for the Fourier transforms, we may conclude that G + can be continued analytically 
into the upper complex uj half-plane and G_ into the lower half-plane. For uj in the 
domain of analyticity of G + or G_ one has the Cauchy integrals, 

n , x 1 f + °° G + (x u x 2 ,uj') , 

G + (x 1 ,x 2 ,w) = — — —duo , 

2m J-oo uj' — uj — lU 

G-(x 1 ,x 2 ,uj ) — / dui . 2.50 

Using the identities (2.47) and (2.50) we obtain the Hilbert transforms 

T , P f+°° ReG±(x u x 2 ,uj , ) J , 
lmG ± (xi,x 2 , uj) = =F— / aw, 

7T J-oo UJ — UJ 

T5 n ( \ ± P f + °° I™G±(x 1 ,x 2 ,uj') , 
ReG±(xi,x 2 ,uj) = ±— = duj . (2.51) 

7T J-oo UJ' — UJ 

In analogy with similar equations in Optics, (2.51) are called the dispersion relations 
for the retarded and advanced Green's functions. 

The Fourier transform of the causal Green's function has the asymptotic form 

G(xi, x 2 , uj) = — S(xi — x 2 ), (uj — > ±oo). 

UJ 

The propagator is a linear combination (2.40) of G + and G_, but since these have 
non-intersecting domains of analyticity for complex uj there are no simple dispersion 
relations for the Fourier transform of the causal Green's function. Even so, there are 
dispersion relations for its diagonal element 

G(x,uj):=G(x,x,uj). (2.52) 

To derive these we first note that it follows from (2.10) that 

C_ (x,x,0) > 0. (2.53) 

From this and (2.31) we obtain the inequality 

/+oo 
K + (x,x,uj)e~^duo > 0, 
-oo 

valid for all (3 G M + . Hence, 

K(x, uj) := K + (x, x, uj) > 0. (2.54) 
From the spectral representation (2.41), by means of (2.47), we obtain 

ReG(^) = gr ImG ^'> f^W, 
V ' ; 7T J-oo uj 1 -uj \e? ±1) 
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lmG(x, uj) = ~K(x, uj) (l ± e' 13 ") . (2.55) 

These are the dispersion relations for the causal Green's function or the propagator. 
Rewriting (2.55) separately for the Bose and Fermi cases, we obtain 

ReG(x, uj) = — J u> u tanh 1 2 J ' ( Bose ) 

ReG(x, uj) = — J — - coth ( -j- J dcu'. (Fermi) 

The equality (2.55) gives 

2lmG(x,u) 

K M = ~ 1±e - P . ■ 

Using this, we may write the density of particles for a system in equilibrium as 

1 f+°° lmG(x,u) , 

^—nL-^ir**- (2 - 56) 

Recall that for equilibrium this density does not vary with time. 

Another form of the dispersion relations for the causal Green's function can be 
expressed by means of the diagonal element of the spectral function (2.44) which we 
denote by 

J (x, uj) :— J (x, x, uj) , (2-57) 

which, as follows from (2.54), is a real valued function. Then, from the spectral repre- 
sentation (2.46), we obtain the dispersion relations 



r 

ReG{x,u) = ~P 



+°° J(x, u') duj' 

Uj' — UJ 2"7T ' 



lmG(x, uj) = --[l± 2n(uj)]J(x, uj), (2.58) 
2 

where the equality 

1 ± 2n(uj) = —, 

v > e ^ =F 1 

has been used. 

If we make use of (2.57), the expression (2.56) for the particle density takes the 
simpler form 

r +oc duj 
p(x,t)= / J(x,uj)n(uj)— . (2.59) 

J-oo An 

Other observables can similarly be expressed by means of integrals involving the spec- 
tral function (2.57). 

For equilibrium systems there is still one more very important type of Green's func- 
tion - the temperature Green's functions, G tem , the definition of which is analogous to 
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that for the causal Green's function except that imaginary time variables are employed. 
Set 

t j = ~ iT j^ U = 1 , 2 ,--) 

with Tj in the interval < Tj < f3. The definition of the one -particle temperature 
Green's function then takes the form 

G tem {\2) := -i [Q(t)C + ( Xi ,x 2 , -it) ± e(-r)C_(a; 1 , x 2 , -it)} , (2.60) 

where r := T\ — r 2 = i{t\ — t 2 ). 

The equilibrium condition (2.30) for the temperture Green's function then becomes: 

G tem (12) | T=0 = ±G tem (12) T=f3 . (2.61) 

Because of this periodicity condition, G tem may be expressed as a Fourier series. Spec- 
tral representations and dispersion relations for the temperature Green's function can 
be derived which are strictly analogous to those for the real-time propagator. 

Which type of Green's function one wishes to use is largely a matter of taste. 
However, in our view, the real-time causal Green's functions are the most convenient 
and so, normally, we prefer to deal with them. 

2.4 Evolution Equations 

The previous two sections were devoted to the study of Green's functions for systems of 
identical particles in equilibrium. We now consider arbitrary statistical systems which 
may or may not be in equilibrium and obtain the equations of motion for the various 
Green's functions. These equations are practically the same for Advanced, Retarded 
and Causal Green's functions. For the sake of concreteness we shall work mostly with 
Causal Green's functions that is with propagators. 

Assume that the system is specified by the Hamiltonian: 

H{t) = J ^(x, t) [K(V) + U(x, t) - /j] tj}{x, t)dx+ 

+- J ip\x, t)ip\x', t)V(x, x')ip(x', t)ip(x, t)dxdx', 

in which K(V) is a kinetic-energy operator. Introduce the effective chemical potential 

/i(l) := fi-Uix^t,), (2.62) 

which incorporates the external fields in U(l), and define the interaction potential 

$(12) :=^(x 1 ,x 2 )6(t 1 -t 2 ), (2.63) 

where $(x 1 ,x 2 ) ■= \[V(x 1 ,x 2 ) + V(x 2 ,x 1 )]. 

In what follows we shall occasionally employ propagators with coincident variables. 
In order to avoid ambiguity, we generalize the definition of the chronological operator in 
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(2.1) to the case of two field operators with the same argument and adopt the following 
convention: 



f ^(1)^(1) := ±^(1)^(1) f ^(1)^(1) := ^(1)^(1). 



(2.64) 



A product of two field operators in which the creation operator is to the left of the 
annihilator has been called the normal form. Thus, the chronological operator by 
acting on the product of ^(1) and reduces it to normal form. This language 

occurred first in connection with the well-known theorem of Wick . 

The equation of motion for the one-particle propagators follows from the Heisenberg 
evolution equation for the field operator and the equality 



d_ 

It 



(±t) = ±5 (t) 



This implies 



d 



G(12) T i j $(13)G 2 (1332)d(3) = 5(12), (2.65) 



where 



5(12) = 5(x 1 -x 2 )5(t 1 -t 2 ). 

A quantity which plays an important role in the theory of Green's functions is the 
self-energy, £(12), which is defined implicitly by the equation 



j E(13)G(32)d(3) =±J $(13)G 2 (1332)d(3). 
Then (2.65) can be transformed to the form 



(2.66) 



G(12) - J S(13)G(32)d(3) = 5(12). (2.67) 



Defining the inverse propagator 



G-\12) := 
we may write (2.67) as 



5(12) -E(12), 



G- 1 (13)G(32)rf(3) = 5(12). 
It then follows from ( 2.66) that 



£(12) =±if $(13)G 2 (1334)G'- 1 (42)rf(34). 



(2.68) 



(2.69) 



(2.70) 



Suppose that G (12) is any propagator for which there is an inverse, G 1 (12), 
satisfying the equation 

j G 1 (13)G (32)rf(3) = 5(12), (2.71) 
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then (2.69) can be transformed into the integral equation 



G(12) = G (12) + J S(13)G(32)d(3), (2.72) 

where the kernel 

5(12) := J G (13)[G' 1 (32) - G _1 (32)]d(3). (2.73) 

The equation (2.72) is called the Dyson equation which is thus a transform, via (2.69), 
of the equation of motion (2.65). 
If we choose 



G \l2) : 

then the kernel (2.73) becomes 

.'(12) = J G (13)[S(32) -S (32)]d(3). 



S(12)-£ (12), (2.74) 



B( 



Since (2.70) and therefore -8(12) involve G 2 , Dyson's equation does not by itself 
determine the one-particle Green's function (2(12). The equation of motion for G 2 
involves G3 and so on up to Gn • There is thus a hierarchy of equations for Green's 
functions as there is for Reduced Density Matrices. In order to short-circuit this chain 
of equations one makes an ansatz for the Green's function of some chosen order thus 
introducing approximations. 

Perhaps the most common such an ansatz is to assume that 

G 2 (1234) = <Sf F (1234) := G(14)G(23) ± G(13)G(24), (2.75) 

which is called the Hartree-Fock approximation. Making this ansatz will lead us to 
equations determing the 1-propagator. 

The Hartree-Fock ansatz enables us to obtain the following expression for the self- 
energy (5.70), 

Ehf(12) = ±i<J(12) J $(3)G(33)d(3) +i$(12)G(12), (2.76) 

where 

G(ll) := lim lim G(12), (2.77) 
in accord with the convention (2.64). Introducing the Hartree potential 

Vff(l) := / $ (12)p(2)d(2), (2.78) 

in which 

p(l) =±iG(ll) =p(x 1 ,h) 
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is the density of particles, we may reduce (2.76) to 



Z HF (12) = V H (1) S (12) + i$ (12) G (12) . 



There is a more general ansatz which appeared in connection with the theory of 
superconductivity and is now usually referred to as the Hartree-Fock-Bogolubov ap- 
proximation. This implies that the two-particle propagator decouples as 



where F and F + are the anomalous Green's functions defined in (2.15). This approxi- 
mation is used not only for superconductors but also often for discussing Bose systems 
with condensate. Under assumption (2.79) the self-energy (2.70) takes the form 



In addition, we will need equations of motion for the anomalous Green's functions which 
would enable us to obtain an approximation to the one-particle Green's function. 

In these two examples we made an assumption about the nature of the second-order 
propagator which enables us to obtain an approximation for the first-order propagator. 
In order to approximate the two-particle propagator one must decouple it from the 
third and higher order propagators by an appropriate ansatz. We now derive the Bethe- 
Salpeter equation and the Brueckner approximation which are two widely used methods 
of moving beyond the 1-propagator and the Hartree-Fock or Hartree-Fock-Bogolubov 
approximations. 

The equation of motion for G2, which involves G3, is as follows: 



Gf FB (1234) = G(14)G(23) ± G(13)G(24) =f F(12)F+(34), 



(2.79) 




. d 



AT(Vi) +//(!) G 2 (1234)^ / $(15)G 3 (125534)rf(5) 



= <S(14)G(23) ± <S(13)G(24). 
Using (2.68) the first term in this equation can be rewritten as 



(2.81) 



,d_ 



tf(Vi) +//(!) G 2 (1234) 




Operating on (2.81) by G, leads to the integral equation 



G 2 (1234) = G(14)G(23)±G(13)G(24)+i / G , (15)$(56)G , 3 (256634H56) 




(2.82) 
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Unhappily this involves G3, otherwise it would be an equation determining G 2 . The 
standard procedure is to "decouple" by an ansatz for G 3 in terms of G — Gi and G 2 . 
In analogy with (2.75), we assume that 

G3(123456) = G(16)G 2 (2345) ± G(15)G 2 (2346) + G(14)G 2 (2356), (2.83) 

which is an expression with appropriate symmetry properties. Substituting (2.83) in 
(2.82) leads to the integral equation 

G 2 (1234) = G 2 , (1234) +i J G(15)$(56)G'(26)G' 2 (5634)d(56)± 

± J G(15)E(56) [G 2 i (2634) - G 2 (2634)] d(56). (2.84) 

In this we used the definition (2.66) for the self-energy and the notation 

G°(1234) := G(14)G(23) ± G{13)G{24). (2.85) 

Since the difference between G 2 and G® arises from interaction between particles, the 
last term in (2.84) includes the interaction potential $ to the second degree and higher. 
Dropping the third term on the right-hand side of (2.84) leads to the Bethe-Salpeter 
equation 

G 2 (1234) = G° 2 (1234) + i J G(15)$ (56) G (26) G 2 (5634) d (56) , (2.86) 

which includes interaction in first order. 

One could try to solve this equation directly or adopt an approach due to Brueckner 
by introducing the so-called T-matrix which is defined implicitly by the relation 

$ (12) G 2 (1234) = Jt (1256) G° 2 (5634) d (56) . (2.87) 

Then from (2.86) and (2.87), using the property 

/ G° 2 (1256) G' 1 (53) CT 1 (64) d (56) = 5 (14) 5 (23) ± 5 (13) 5 (24) , 

it is possible to show that 

T (1234) = $ (12) 5 (13) 5 (24) + i$ (12) J G (15) G (26) T (5634) d (56) . (2.88) 

The approach based on (2.88) is called the Brueckner approximation. If (2.88) can be 
solved for T, then (2.87) can be used to obtain G 2 . 

In concluding this section dealing with the equations of motion, we shall prove 
an important - even somewhat startling - result which echos the famous theorem of 
Kohn-Hohenberg. From (2.65) assuming that t 1 7^ t 2 , we have 



J $(13)G 2 (1332)d(3) 



iA_ X(Vl ) +//(!) 



G(12). 
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If, as in the standard case, the interaction betwen particles is symmetric so V(x±, x^) = 
$(x 1 ,x 2 ), then the average of the Hamiltonian (2.25) can be cast in the form 



< H >= ± 



- / lim 

2 J (21) 



iA + ^( Vl )-^(i) 



G(12) dx x . 



^2.89) 



Further, the internal energy (2.26) becomes 



E = ±- / lim 

(21) 



d 

i— + K(V 1 ) + U(x u t 1 )+pL 



G(12)d Xl . 



(2.90) 



For equilibrium, when U(x,t) = U(x), we may invoke the Fourier transform in (2.90) 
to obtain 



E = ±- 



//" + 00 fji^ 
dx [lu + K(W) + U(x) + fj] e +tuJ °G(x, u) — , (2.91) 



where the notation (2.52) is used. Finally, using the spectral representation (2.46) we 
obtain from (2.91) 

E = - dx [cu + K(V) + U(x)+ii] J(x,u)n(cu) — , (2.92) 
2 J J-oo 2tt 

where J(x,u) is the diagonal element (2.57) of the spectral function (2.44). 

This means that in order to calculate the internal energy of a system in equilibrium 
we need to know only the 1-propagator, or the diagonal of the corresponding spectral 
function. Note that the expressions (2.91) and (2.92) are exact. Of course, it will only 
be for moderately simple models that we shall be able to obtain the exact 1-propagator. 



2.5 Wave Functions 

In this section we shall introduce a family of 1 -particle functions which can be thought 
of as Fourier transforms of time-dependent natural orbitals. Since, in quantum me- 
chanics, time and energy are conjugate observables the effective wave functions which 
we shall define will be associated with 1-particle energies. The effective Hamiltonian 
of which these wave functions are eigen-functions contains the self-energy operator and 
therefore, in analogy with Hartree-Fock orbitals, they may be thought of as describing 
"dressed" quasi-particles which encapsulate a major portion of the correlation between 
the bare particles. As far as we are aware the precise relation between natural orbitals 
and the effective wave functions has never been studied so we do not assert that the 
effective wave functions in the limit specified in (2.9) become the natural orbitals of 
the Hartree-Fock approximation to the state of the system. We do not rule this out as 
a possibility but will be surprised if it proves to be the case since, while the Hartree- 
Fock approximation and the method of effective Hamiltonian of this section are both 
attempts to take into account a portion of the interaction energy, they are probably 
not strictly equivalent. We believe that this is an issue which deserves further study. 
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For a system in equilibrium we reformulate the equations of motion in terms of 
their Fourier transforms. Thus the inverse propagator can be expanded as 



G-\ Xl ,x 2 ,u)e-^ , (2.93) 

-oo Z7T 



in which as in (2.39) we have used the notation t = t\ 2 . Then the evolution equation 
(2.69) takes the form 

J G^ 1 (xi,x 3 , uj) G (x 3 , x 2 , uj) dx 3 — S(xi — x 2 ) . (2.94) 

Similarly for the self-energy, 

/+°° duj 
E( Xl ,x 2 ,u) — . (2.95) 

Introducing the notation 

H(x, x, uj) := [K(V) + U(x)]5(x - x ) + E(ar, x', uj), (2.96) 

for the Fourier transform of the inverse propagator (2.68) we find 

G-\x, x', uj) = (uj + fi)5(x - x') - H(x, x', uj). (2.97) 

With respect to a complete ortho normal set of one-particle functions, {d n (x)}, we 
have expansions 

G{x,x',Uj) = Y,G m n(u>)# m (x)r n (x') , 
m,n 

G-\x, x', w ) = £ G^ n (uj)# m (x)r n (x>) , (2.98) 

m,n 

and hence, by (2.97), 

G^ n (uj) = (uj + fi) 5 mn - H mn (uj) , (2.99) 

where 

H mn (uj): = f $* m (x)H(x,x',uj)$ n (x')dxdx'. (2.100) 
The equation of motion (2.94) becomes 

Y,G m 1 p (oj)G pn (u) = 5 mn . (2.101) 

p 

Similarly, we can express the spectral function (2.49) in the form 

J(x,x',Uj) = J2 J mn(0j)^ m (x)^ n (x'), (2.102) 
m,n 

with 

Jmn (w) = i [G mn (w + i0) - G mn (uj ~ i0)] . 
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In terms of these we may write the corresponding spectral representations and disper- 
sion relations. 

Expressions such as the above are most useful if G mn (uj) is a diagonal matrix. We 
define an effective Hamiltonian, H(x,u), as an integral operator such that 

H(x,uj)<f(x,uj) := J H(x, x',u)ip(x',u)dx', (2.103) 

with kernel (2.96). From (2.94) and (2.97) the equation of motion then takes the form 

[u + fi- H(x, u)]G(x, x', u) = 5(x - x'). (2.104) 

If we substitute (2.98) in this and assume that ~& n can be chosen so that G mn {u) = 
G m (u)5 mn , we find 

J2 G n {u) [u + fi- H(x, u)]& n (x)&* n (x') = S(x -x')=J2 Mx)K(x')- 

n n 

In the last step we have used the assumption that is a complete orthonormal set. 
This equality will be possible only if 

G n (uj) [u + /I - H (x, uj)] i?„ (x) = § n (x) . 

In particular it is clear that the $ n are eigenfunctions of H(x, u>) with eigenvalues 
depending only on uj and the constant \i. Set 

H{ X) u)d n {x) = H n {uo)d n {x). (2.105) 

Equation (2.105) is called the effective wave equation and the functions i? n , effective 
wave functions. Some authors also call (2.105) the effective Schrodinger equation. 
Comparing (2.104) and (2.105), we see that the poles of the propagator are defined by 
the eigenvalues of the effective Hamiltonian. The eigenfunctions, q 9 n {x) ) describe quan- 
tum states, indexed by n, analogous to the familiar wave functions for the 1-particle 
Schrodinger equation. However, the effective Schrodinger equation is an integral rather 
than a differential equation and contains the additional parameter oj. 
Because of (2.105), the matrix elements (2.100) take the form 

Hmn M = H n (uj) 5 mn . (2.106) 

The Fourier coefficients of the inverse propagator (2.99) then become 

G^ n (Lu) = [cu + fi- H n (lu)} 5 mn . (2.107) 

Substituting (2.107) into the equation of motion (2.101), we obtain 

G mn (to) = G n (cu)5 mn , G n (u>) = , ■ . (2.108) 
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Equation (2.108) defines G n (uj) as a function of uj on the whole of the complex plane 
apart from the real axis. For real values of uj, G n is defined by means of the spectral 
representation (2.46). This, together with the spectral function J mn (uj) = J n (uj)5 mn , 
yields 



/+oo 
J n (u) 
-oo 



1 ± n(uj') 



n(uj') 



uj — ui' + iO uj — uj' — iO 
When the spectral function, J n (ui), has no poles on the real axis, (2.109) implies that 

J n (uj) = i[G n (uj + iO) - G n (uj - iO)]. (2.110) 

Indeed, in the expression in the square bracket of (2.109), the terms involving n(uj') 
cancel and leave an integrand J n (uj)/(ui — uj') to which Cauchy's formula applies. 

In order to calculate J n (u), we need to know how (2.106) depends on uj, so in view 
of (2.96) we separate the self-energy into two parts, which are respectively independent 
and dependent on uj, by defining A(x,x',uj) to satisfy 



E(x, x', uj) = E(x, x' , 0) + A(x, x' , uj) 
Then , in view of (2.103) and (2.105), 

H n (uj) = H n (0)+A n (uj), 
where 

H n (0) = J r n (x) [K(V) + U(x)} d n (x)dx + J d* n (x)^(x,x', $)d n (x')dxdx' , 

and 



(2.111) 



(2.112) 



A n (uj) = J $* n (x)A(x,x',uj)$ n (x')dxdx'. 



The second term in (2.112) can be related to the propagator G n (uj) by means of the 
equation of motion and has analytical properties similar to those of G n (uj). Just as 
J n (uj) was defined in connection with the spectral representation of G n (uj), so we may 
introduce a spectral function, T h (uj) such that for complex values of uj, with lm(uj) ^ 0, 



A 



n(t0) = J 



■+°° T n (uJ r ) duj' 
-oo uj — uj 1 2n 

For real uj, the spectral representation is 

1 ± n(uj') n(u>' 



(2.113) 



/+0O 
T n (uj') 
-oo 



uj — uj' + iO uj — uj' — iO 
In analogy with (2.110), it can be proved that 



i [A n (uj + iO) - A n (uj - iO)] 



(2.115) 
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Thus, when we take (2.112) and (2.113) into account, the propagator (2.108) as- 
sumes, for complex uj off the real axis, the value 

G„(u) = L - ff.(0) + y. - r M^}~ ' . (2 .n6) 

[ J-oo UJ — UJ An ) 

In preparation for applications, we define an effective energy as 

E,H:=H,(») + W W ^^ < 2 ' 117 > 

J-oo UJ — UJ 271 

Then substituting (2.116) in (2.110) we obtain the following important and very useful 
expression for the spectral function 

Uu) = (2 - 118) 



The expression 

J-oo UJ — UJ 271 



ReA 

is called the energy shift, and 



r » M = -r±^U (2 ' 120) 

is called the decay rate for the state ~& n {x). 

Neglecting the decay is equivalent to retaining only the real part of the poles of the 
propagator. In this case, (2.110) together with (2.47) imply that 

J n (u) = 2n5 {uj - E n (u) + n) . (2.121) 

When a smooth function f{ui) has only simple zeros, u>i, that is f'iuji) ^ for all i, 

xfff \\ \- S(u - Uj) 

,(/(u)) = ?TOT- 

Applying this formula to (2.121) yields 

where uj ni are defined by the equation 

uj ni = E n (uj ni ) - fi, (i = 1, 2, . . .) 

and 

KM = ^*M; KW ^ i- 
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If we use the spectral function (2.122), the propagator (2.109) takes the form 



G, 



1 ± n{u ni ) n(uj ni ) 



UJ 



0J ni + i0 UJ - UJ, 



i0 



(2.123) 



which is simpler and easier to use than when the spectral function (2.118) is inserted 
in (2.109). But one should remember that this last form is appropriate only if decay 
processes may validly be neglected. 



2.6 Scattering Matrix 

In his very first paper about Quantum Mechanics, Heisenberg presented a theory in 
terms of infinite matrices which, as far as possible, involved only observable properties 
of the physical situation under examination. Later, he introduced his scattering matrix 
appropriate for a situation in which one or more particles came into a "black box" and 
emerged in a different direction or as different species. He defined a matrix which con- 
nected the intial to the final state of the system. In this section, we use propagators to 
take into account the nature of the "black box" and to calculate the scattering matrix. 
Our Soo, defined in (2.139), is a particular case of Heisenberg's concept. The ideas 
of this section play an important role in providing a practical method of calculating 
Green's functions. 

We envisage a system which is evolving under the control of a Hamiltonian H and 
ask how it behaves if the Hamiltonian is modified by the addition of a term, AH, which 
may or may not be "small". Both H and AH are functionals of the field operators, 
which we indicate by 

H = H{^}, AH = AH {4)}. (2.124) 

The new Hamiltonian is 

H A := H + AH, (2.125) 

where 

H = H{4> A }, AH = AH{4> A } (2.126) 

are the same functionals as in (2.124) but functionals of new field operators ip A which 
satisfy the Heisenberg equation 

d 

(x,t) = [ip A {x,t),H A (t)\ 
Thus, the time dependence of the new field operators is given by 

V'a {x, t) = U+ (t) iP A Or, 0) U A (t) , (2.128) 
where the evolution operator satisfies the equation 

i^-U A (t) = H A U A (t) , 



(2.127) 
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with H A and H A (t) related, analogously to (2.128), by 

H A (t) = U A -{t)H A U(t). (2.129) 

Since the statistics of particles is independent of the Hamiltonian, the new field oper- 
ators satisfy the same commutation relations as the old ones. In other words, 

[i/j A (x,t),ijj A (x',t)] lz = 0, 

and 

iPa(x, t),ip A (x', t) _^ = 5(x-x'). 

Since H will be assumed to be self-adjoint and we want the creation operators 
4> A (x, t) to satisfy the same equation (2.127) as ip A (x, t), we shall also assume that AH 
is self-adjoint. Thus 

AH + = AH. (2.130) 
This imposes a restriction on AH so that if, for example, 

AH :=- J ip A (x)Api(x, x', t)tl) A (x')dxdx' , (2.131) 

in which ip A (x) := ip A (x, 0) and Ajj, is a complex-valued function, then (2.130) implies 
that 

Afj,*(x, x', t) = A/i(x, x', t). (2.132) 

Corresponding to the new field operators, new propagators can be defined just as 
the old ones were in (2.3), (2.11), and (2.13). Thus, for the one- and two-particle 
propagators we have, respectively 

(5(12) := -i < TVa(1)V4(2) >, 

G 2 (1234) := - < ^(1)^(2)^(3)^(4) > . (2.133) 

The statistical averaging in (2.133) does not presuppose that the system is in equil- 
brium. It is taken with respect to the system density matrix, pit), prescribed by p(0) 
at initial time zero. 

As compared to the equations of motion for the old propagators, those for the new 
ones, (2.133), contain an additional term which results from the commutation of ip A 
with 

AH (t) = Ul (t) AHU A (t) . 
For the one-particle propagator G, the equation of motion takes the form 

i^-G (12) = 5 (12) i <L T [ip A (1) , H A (t,)} ^ (2) > . (2.134) 

Oti 

If we take the additive term in the form (2.131) we find 

ip A (x, t) , AH (t)] = - J Afi(x, x', t)ip A (x', t)dx', 
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and obtain 



iA_/r( Vl ) +1,(1) 



G (12) + J Afi ( Xl ,x 3 , h) 5 (h - t 3 ) G (32) d (3) =f 

=F i J $ (13) G 2 (1332) d(3) = 6 (12) . (2.135) 
Notice that when A fx — > 0, so Aff tends to zero, 

■0a (x, -> V (a:, , H A ^H -> if, 

while (2.135) becomes (2.65). Since we have not assumed that Aii is small, the equa- 
tion (2.135) for the propagator is valid for any AH of the form (2.131). 
We now introduce the scattering matrix 

S(t 1 ,t 2 ) := f exp j-i jf* 1 Afi (t) dij , (2.136) 

in which 

Afi (t) = U + (t) AHU (t) , 

and the term 

AH :=- J ^ (x) A/i (x, x, t) i> (x') dxdx (2.137) 

has the same form as AH in (2.131) but contains the original field operators tp and ^ . 
The scattering matrix, S(ti,t 2 ), has the following properties: 

S(t,t) = l, lim S(t u t 2 ) = l, 

lAti — ►U 

S(t 1 ,t 2 )S(t 2 ,t 3 )=S(t 1 ,t 3 ), (2.138) 

s+ (h,t 2 ) = S- 1 (t 1 ,t 2 ) = S(t 2 ,t 1 ). 

The limits of integration, t\ and t 2 in (2.136), can be any real numbers. In particular 
we may define the balanced double limit 

: = lim S (t, -t) = S (+00, -00) , (2.139) 

giving 

Sec := f exp j-i y °° Aif (i) dt} . (2.140) 

By considering the time-ordered product of with the field operator VK-Qj we infer 

TS^ (1) = 5 (+00, h) i> (1) 5 (t u -00) . 

In a similar manner, we can obtain the time-ordered product of with several 
field operators. For example, 

f 5^(1)^(2) = 0(t 12 )5(+cx),t 1 )V(l) 1 S(t 1 ,t 2 )^ t (2)5(t 2 ,-oo)± 
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±e(t 21 )S(+oo, t 2 )^(2)S(t 2 , t^i^Sih, -oo). 

With the help of the scattering matrix we can express the new propagators in 
(2.133) by means of the old field operators: 

~ ; <T5^(l)^t (2) > 
G (12) = -i , 

G 2 (1234) = - < ^W^ ,(3 ^ t(4 >. (2.141) 

Similarly, the propagator, G n , of order n, can be expressed by means of the field 
operators as 

a. (i . . . 2») = (_r <^l')^Mf("H)..i'M> (2 142) 

When *0 in (2.13) is replaced by ip A , we obtain (2.142). The proof is obtained by 
considering their equations of motion. We illustrate this in the case of the first order 
propagator in (2.141). Differentiate G with respect to t\ and take into account the 
relations 

i—S (+oo, h) = -S (+oo, tO AH (*0 , 
oti 

d 

i—S(t u t 2 ) = AH (tO S (ti,t 2 ) , 

C/ti 
Oti 

i^-S (t x , -oo) = AH (tO 5 (fx, -oo) . 

Oti 

From (2.141) we then find 

i^-G (12) = f (12) - fS ™ & (1) • " + AH »' (2) > . 

dti < Soo > 

This is seen to reduce to (2.135), when we recall the definition of H(t) and note that 
(2.137) implies that 

(x, t) , AH (t)] = - J Afi(x, x', t)ip(x, t)dx. 

It is important to realize that the representation (2.141) for the propagators (2.133) 
is valid for arbitary shifts AH of the Hamiltonian (2.125) whether or not the system is 
in equilibrium. It is frequently assumed that H is the Hamiltonian for non-interacting 
particles while AH describes their interaction. In this case, (2.141) is referred to as 
the interaction representation for the Green's functions. 
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Consider a particular case for which AH is so small that we may restrict ourselves 
to terms which are linear in 

AG (12) :=G(12)-G(12). (2.143) 
For such an approximation, linear in AH, the scattering matrix (2.140) is 

/+oo 
AH(t)dt. (2.144) 
-oo 

From (2.141), it follows that 

r+oo 

G(12) ~G(12)+iG(12) / <TAH(t)>dt- 

r+oo A 

- / < TAH (t) i) (1) V f (2) > dt. 
If AH has the form (2.137), then 

< TAH (t) >= =f J A/j, (x 3 , x 4 , t) 5 (t 3 -t)8 (t 4 — t)G (43) d (34) , 

< TAH (t) i) (1) ^ (2) >= ± J A/i (x 3 , x 4 , 5 (t 3 - t) 5 {U - t) G 2 (1432) d (34) . 
Therefore the difference (2.143) becomes 

AG (12) = ± J [G (12) G (43) - G 2 (1432)] Afi (x 3 , x 4 , t 3 ) 5 (t 3 - f 4 ) d (34) . 

This last expression, with the notation 

Afi (12) := Afi (xi, x 2 , h) 5 (h - t 2 ) , (2.145) 

takes the form 

AG (12) ~ ± J [G (12) G (43) - G 2 (1432)] Afj, (34) d (34) . (2.146) 
Define the variational derivative 

- = lim AG ( 12 ) (2 U7) 

5/i (34) ' a^o A/i (34) ' v • ; 

Then (2.146) yields the variational representation 

G 2 (1234) = G (14) G (23) =F (2.148) 
for the two-particle propagator. 
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If A/j, is diagonal in x, so that 

A/x (x, x', t) = A/j, (x, t) 5 (x — x') , 

then (2.146) gives 

AG (12) ~ ± J [G(12)G (33) - G 2 (1332)] Afi (3) d (3) , (2.149) 

where 

A/j(l) = A/j (a;i,ti). 
Defining the variational derivative 

SG (12) , AG (12) 
= , / := lim . \A 2.150 
5/i(3) a^o A/i(3) ' v ; 

we obtain the variational representation 

G 2 (1332) = G (12) G (33) =f (2.151) 

for the two-particle propagator. Note that this has precisely the form which we need 
in the evolution equation (2.65). 

2.7 Perturbation Theory 

Using the definition of the scattering matrix it is straightforward to develop perturba- 
tion theory for Green's functions. 

Suppose that the unperturbed Hamiltonian in (2.124) is H = H with 

H = J V> f (x) [K (V) - /i (x, £)] i> (x) dx, (2.152) 

where \x (x, t) may include the chemical potential, external fields and some self- 
consistent mean-fields. It might, for example, include the Hartree potential (2.78). 
We assume a perturbation, corresponding to a two particles interaction, in the form 

AH = V f (x) ^ (x') V (x, x') i\) (x') i\) (x) dxdx'. (2.153) 
Denote the unperturbed Green's function by Go- Employing the commutator 

[^(l),//o(ii)] = {^(V 1 )-/i(l)}V(l), 
we obtain the equation of motion for the unperturbed propagator, Go: 

i^--K (Vi) + n (1) G (12) = 5 (12) . (2.154) 
at i 
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This equation can be solved by the method of expansion in wave functons, discussed 
in Section 2.4. 

We proceed as follows in order to obtain perturbative corrections to Gq. Expand 
the scattering matrix (2.140) in the series 



\k „,„ k 



S oo = ^ + J2-TT Tl[AH(t p )dt p . (2.155) 
k=i K - J ~°° P =i 

Truncating this series at order k and substituting it into the interaction representation 
(2.141) gives rise to the k-th order approximation, G^ k \ to the propagator G. This 
approximation includes powers of AH up to the A;-th. Clearly, the zero-th order 
approximation is the unperturbed propagator 

G (0) (12) := Go (12) . 

The k-th order approximation, G^ fc , for the n-particle propagator, G n , can be obtained 
in the same way. 

Taking account of the equality 

f + °° < TAH (t)>dt=- [ G { 2 ] (1221) $ (12) d (12) , 

the first-order approximation for the one-particle propagator becomes 

G (1) (12) = G (12) + 

+ % J [Gf (134432) - G (12) G^ 0) (3443)] $ (34) d (34) . (2.156) 

Here, the n-particle propagators, G^\ satisfy equations of motion involving the un- 
perturbed Hamiltonian, H . 

The equation of motion for G^ is as follows: 



iA_ x(Vl ) + IM (i) 



G { 2 ] (1234) = 5 (14) G (23) ± 5 (13) G (24) . (2.157) 



Since the solution of a first-order differential equation is determined by the intial values, 
by comparing this equation with (2.154) we may conclude that 

Gf (1234) = G (14) G (23) ± G (13) G (24) . (2.158) 

In a similar manner, from the equation of motion for G^ together with (2.154) and 
(2.157), we conclude that 

GP (123456) = G (16) G^ (2345) ± 

± G (15) G { 2 ] (2346) + G (14) Gf (2356) . (2.159) 
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From (2.158) it follows that G^ can be expressed as a sum of six terms involving only 
the first-order propagator G . We symbolize this by 

Gf = £(±l) p G G Go, 



where P runs over the group of permutations of {4,5,6}. By (±l) p we mean +1 for 
bosons and, for fermions, the parity, +1 or —1, of the permutation P. Each G must 
contain a variable from the first three and from the second three variables in The 
terms which occur can be obtained from Go(16)Go(25)Go(34), in which {1, 2, 3} are in 
natural order and {4, 5, 6} are in reversed order, by the six possible permutations of 
{4,5,6}. 

Similarly, any n-particle unperturbed propagator can be presented as a sum 



Gi°) = ]T(±l) P GoGo...G c 



(2.160) 



The sum will contain n\ terms obtained by a rule which is the obvious generalization 
of the preceding example. In quantum field theory, this result is usually referred to as 
Wick's theorem. 

Using (2.158) and (1.159) in (2.156) leads us to 

G (1) (12) = G (12) ± % J G (13) $ (34) Gf (3442) d (34) , (2.161) 

for the first-order approximation to G(12). The second-order approximation is given 
by 

1 r+oo 



G (2) (12) = G (1) (12) + - / 

2 J— oo 



f {[G (12)+ #(1)^(2) 
x [AH(t 3 ) - 2 < AH(t 3 ) >] AH(t 4 )} > dt 3 dt 4 , 



x 



in which we need to substitute (2.153) and make use of (2.160). 

In a similar manner we can obtain any fc-th order approximation for an n-particle 
propagator G!£\ 

Another approach to perturbation theory uses the Dyson equation (2.72). Setting 



S (12) = 0, 
we find that the inverse propagator (2.74) is 



Go 1 (12) 

and the kernel (2.73) becomes 



iJL-KW+viQ.) 



5(12), 



(2.162) 



(2.163) 



£(12) = j G (13) E (32) d (3) 
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The Dyson equation, (2.72), takes the form 

G (12) = G (12) + J G Q (13) E (34) G (42) d (34) . (2.164) 

Successive approximations to G are obtained by iterating (2.164) according to the 
scheme 

G (k) _^ E ( fc +i) _^ G (fc+i)_ ( 2165 ) 
The first-order approximation gives 

G (1) (12) = G (12) + y G (13) E {1) (34) G (42) d (34) , (2.166) 

in which 

E« (12) = ±i J $ (13) Gf (1334) G^ 1 (42) d (34) . 
With G^ given by (2.158), we have 

(12) = i$ (12) G (12) ± i(J (12) y $ (13) G (33) d (3) . (2.167) 

Therefore (2.166) coincides with (2.161). 

Note that (2.167) has the same apparent structure as the Hartree-Fock self-energy 
(2.76) but is quite different since E^ 1 -* involves Go, while Tihf involves G. 

By (2.66) the Dyson equation (2.164) can be transformed into 

G (12) = G (12) ± i J G (13) $ (34) G 2 (3442) d (34) . (2.168) 

There is an equivalent approach to obtaining successive approximations to G in 
which we adopt the recursive scheme 

G (k) ^ G (k) ^ G (fc+i)_ (2.169) 

Substituting G^ from (2.158) in the right-hand side of equation (2.168), again gives 
us (2.161). In order to find G^ 2 \ we need This can be obtained from either (2.141) 
or (2.65), giving 

G { 2 ] (1234) = G i 2 ) ( 1234) + 

/+oo , _ . . , 

< T Gf ] (1234) + rj, (1) ^ (2) V f (3) ^ f (4) AiJ (i) > dt 
-oo 

With the help of (2.153), this yields 

G { 2 ] (1234) = Gf (1234) + 

+ i J [Gf (12566534) - G< 0) (1234) Gf (5665)] $ (56) d (56) . (2.170) 
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The equation of motion for Gf\ together with the Wick's theorem (2.160), yields 

Gf ] (12345678) = G (18) G$ 0) (234567) ± G (17) (234568) + 

+G (16) Gf ) (234578) ± G (15) (234678) . 
Invoking (2.159), we obtain 

Gf ] (12345678) = 
= G { 2 ] (1278) Gf (3456) ± Gf (1378) G { ? (2456) + 
+G ( 2 ] (1478) (2356) + G { ° ] (1256) (3478) ± 
±Gf (1356) (2478) + Gf (1456) (2378) , 

where is given by (2.158). 

The two approaches to perturbation theory, are equivalent to one another, and 
applicable whether or not the system is in equilibrium, as long as AH can be considered 
small. Since Green's functions may be interpreted as generalized functions, convergence 
of the perturbation sequence is to be understood as convergence of functionals, - for 
example, as the convergence of the expected values for some local observables. 

It may happen that terms in the perturbation expansion diverge when they involve 
integrals of products of several propagators. We therefore illustrate, by means of an 
example, two methods, which we call implicit and explicit, of avoiding such divergences 



Consider the Fourier transform of (2.163) for a system in equilibrium, 

Gq 1 (x, x', uj) = [u + n - H (x)} 6{x- x') , (2.171) 

where 

H{x) — K (V) + U (x) (2.172) 
is hermitian. The propagator corresponding to (2.171) can be expanded in the form 

Go (x,x',u) = Y J G n {oo)^n{x)r n {x l ) (2.173) 

n 

with respect to one-particle wave functions satisfying the equation 

H(x)il> n (x) = E n il> n {x), (2.174) 

with expansion coefficients 

GnM= l — (2.175) 
defined on the complex c^-plane. As was explained in Section 2.4, on the real cu-axis 

G.(«0= ^"^U ^ . (2.176) 

nK ' w- E n + i0 W- E n -i0 v ; 



The eigenvalues E n are real since we assumed that H (x) is hermitian. 

Starting from the zero order approximation (2.173) for the propagator, in successive 
approximations we meet a variety of integrals with integrands containing products of 
factors such as (2.176). For example: 

/ + OC fjfj 
e + ^G m (u)G n (u) — , 
-oo 

/ + OO fl/.j 
^ e +i »°u;G m (u,)G n (u;)-, 

/ + OO fjrj 
e + ™°G m (u;)Gg(u;)G n (u;)—, 
-oo Z7T 

/ + 00 ft/,) 
^ e^°uG m (u) G e (w) G n ( w ) — . (2.177) 

The integrals in (2.177) are symmetric with respect to permutation of their indices. 

When two or more of the indices in (2.177) coincide, a direct integration of G^(uj), 
or in general of G^(u), would lead to a divergence. This would be due to the fact that 
the propagator G n {uj) must be treated as a generalized function or distribution and 
products of distributions with coincident singualarities are not uniquely determined. 
They require additional definition. 

There are two approaches, so-called implicit and explicit, which avoid the diver- 
gences [59]. 

(i) For the implicit approach, instead of defining G^(u) directly we replace it by 
the limit of a sequence. We first calculate the integrals in (2.177) for distinct indices, 
giving 

Imn A mn -\- A nm , 
^mn E m A mn -\- E n A nm , 
Imln B m £ n -\- Bg nm -\- B nm £ , 
I 'min = Em B m e n + EgBg nm + E n B nm £, 



where 



and 



n (E„ 



B 



min 



- C/ m - C/ n 

n (E m ) [1 ± n (Eg)} [1 ± n {E n )\ ± [1 ± n {E m )\ n (Eg) n (E n ) 



(E m — Eg) (E m — E n ) 

We then define the integrals (2.177), for coincident indices, by the following limiting 
procedure: 

I nn := lim I mn = -f3n (E n ) [l±n (E n )\ , 

m— >n 

4 := lim I' mn = n (E n ) + E n I nn , 
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Intn '■= lim Imln = (2.178) 

m—*n — 



E n I nn — E(I n £ 

l nln '■— i,™ 1 mJLn 



InPn '■ Hill I m 0„ 



m-n™" E n -E t 

I nnn := lim l nln = Z-n (E n ) [1 ± n (E n )\ [1 ± In {E n )\ , 

i^n Z 

Innn • ^n£n Am "I" E n I nnn . 

t—>n 

That is, we let the indices coincide only after integrating but not in the integrands. 
(ii) For the explicit approach, we first define the k-th power of G n as 

Gl M - 1±nM k =F (2.179) 

We then substitute the expression (2.179) directly into the integrands. If under the 
integral, it is multiplied by a function, f(u) say, which is differentiable of order k or 
more and vanishes sufficiently rapidly at infinity, then integration by parts gives 



DO 



Using the relations 



and 



-j-n (lo) = —f3n (lo) [1 ± n (lo)] , 
dto 

d 2 d 

-n (co) = -p [1 ± 2n (cv)] — n (co) , 



duo 2 1 du 

it is straightforward to show that these two approaches to defining the integrals (2.177) 
with coincident indices are equivalent. 

Thus, by using either perturbation method, we obtain a sequence, {G^}, of ap- 
proximations for a propagator. With these we may then find a sequence, {< A >&}, 
of expected values for an element, A, of an algebra of observables. For any particular 
problem, it will be the convergence properties of the sequences < A > k that determine 
convergence of our procedure. 



2.8 Excited States 

It frequently happens that we are more interested in the spectrum of the quantum 
states of an equilibrium system than in its statistical states. As we proved in Section 
2.4, the spectrum is given by the poles of the propagators in the cu-representation and 
can be found from the effective wave equation (2.105). 
The effective Hamiltonian 

H (x, lo) = H (x) + M (x, lo) (2.180) 
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consists of two terms. The first is the Schrodinger operator (2.172) and the second, 
called the mass operator, is an integral operator the kernel of which is the self-energy, 
E: 

M (x, uj) (p (x, uj) :— J E (x, x , uj) (p (x, uj) dx'. (2.181) 
Thus the effective wave equation is 

[H (X) + M (x, UJ)} if n (X, UJ) = H n (UJ) if n (x, uj) . (2.182) 

This equation differs from the usual Schrodinger equation in three respects: 

1. It is an integro-differential equation; 

2. Since, in general, the effective Hamiltonian, (2.180), is not hermitian, the eigen- 
values, H n (uj), are not necessarily real; 

3. The self-energy is normally defined only by means of the perturbation theory 
described in the previous section. 

Therefore in order to solve (2.182) it will be necessary to use perturbation the- 
ory, but (3) implies that the familiar method of solving the Schrodinger equation is 
not applicable since this assumes that the Hamiltonian is known before recursion for 
successive approximations commences. We therefore attempt a solution by setting 

oo 

x,uj) (2.183) 

fc=0 

in which e is a formal expansion parameter introduced to keep track of the successive 
orders of approximation. At the end of the calculation we let e — > 1, hoping that the 
series converges! Each is an integral operator such that 

M {k) 

(x, uj) if (x, uj) := (x, x', uj) ip (x', uj) dx', 

with the kernel, M^ k \x,x' ,uj), defined by the corresponding approximation for the 
self-energy. In the zero approximation, E (12) = 0, we have 

M (0) (x, x', uj) = E (x, x', uj) = 0. (2.184) 

The first approximation is given by 

/+oo 
E (1) (12)e l ^t, (2.185) 

with self-energy (2.167) and t = t 12 . Employing the expansion (2.173) for the zero- 
order propagator G , we find 

i f +OO e^°G (x,x>,uj)p = ±p ( ( ,) (x,x'), 
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where 

pf ] {x,x') :=Y.n{E m )4> m {x)r m {x') (2.186) 

m 

is the 1 -matrix in zero approximation and Ek and tpk are defined by (2.174). Thus 
(2.185) reduces to the expression 

M« (x, x', u) = ±$(x, x')pf ] (x, x 1 ) +5(x-x') J $(x, x")p (0) ( x")dx", (2.187) 
in which 

(x) := pS 0) fox) = £n(£ m ) (x)\ 2 (2.188) 

m 

is the particle density. Note that (2.187) does not depend on u. 
However, in higher approximations, that is for k > 2, the kernel 

/+oo r -, 
£« (12) - E^- 1 ) (12) e iujt dt (2.189) 
-oo 

will, in general, depend on uj. 

The series (2.183) for the mass operator suggests that we should seek a solution of 
the eigenvalue problem (2.182) in the form 

oo 

k=i 



Vn (x, u) = ^ n (x) + £ eVi fc) (x, • (2-190) 

k=l 

The functions ip^ can be expanded in the basis {ip n }, so that 

^ ) (^^)=Ec«Mt( a ;) ) (2.191) 

m 

and, since we assume that the ip m are orthonormal, 

cS&M = (^ mj rf>), fc>l. (2.192) 

There is an additional condition on these coefficients expressing the fact that we want 
ip n (x,uj) in (2.190) to be normalized, 

(<p n ,<p n ) = l. (2.193) 

Substitute the expansions (2.183) and (2.190) into (2.182) and equate like powers 
of e. In zero order, this gives (2.174) and in first order: 

[H(x) - E n ] <pW(x,u) = [hW(u) - MW(x,u)] Mx). (2.194) 
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This can be considered as a nonhomogeneous equation defining ip£\ Using (2.191) we 
may transform (2.194) to 

cS (w) (£ m - E n ) = 5 mn HW ( W ) _ M « ( W ) , (2.195) 

in which 

MffiM := | r m {x)M^(x,x',u)^')dxdx'. (2.196) 

The coefficients c^|(a;), when m — n, are not specified by (2.195) but can be ob- 
tained from the normalization condition (2.193). If we substitute in this condition the 
expansion (2.191) and equate equal powers of e, we find 

<£(<»>) = (*»,<pS ) )=0- (2 - 197) 
It follows that the first order corrections to the eigenvalues and eigenfunctions are 

rf ) (w) = i;cS2iH^m(x), (2.198) 

rrt 

where 

= (2.199) 

In the second order, we find 

[H{x) - E n \ <pS\x,u>) = [HW(u) - MW(x,a;)] o;)+ 

+ [H&>(u) - M^(x,co)} i> n {x). (2.200) 
Together with (2.191) this leads to 

<&(«>) - £») = Wf ] M + cffi M M - 

- E M S M cff M - M^) ( W ) . (2.201) 

i 

From this we obtain 

c mn KW)-Z^, E _e)(E,-E) (E — E ) 2 E -E ' 1 ; 

where in this last equation m ^ n. As before, c^(u;), with coincident indices, cannot 
be obtained from (2.201), but can be deduced from the normalization requirement 
(2.193) which, to within a phase factor, implies 

2 

(2.203) 



m^n 



Ml!2 



E - E 

J -'m J -'n 
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Another method of defining the second-order correction, H^(u), is to multiply 
(2.200) by ipn( x ) an d integrate the resulting equation with respect to x. Noting that 
H(x) is hermitian and using (2.107) leads to 

H<? ( W ) = (>n, M^M) + M£) ( W ) , (2.204) 

where </?^ is defined as the solution of the nonhomogeneous equation (2.194). Thus 
the correction, is expressed by means of quantities with the common label n. At 
first sight, (2.204) looks simpler than (2.202) which involves a summation. However, 
this second approach requires the solution of the differential equation (2.194) in order 
to evaluate the scalar product on the right hand side of (2.204). Which is the optimal 
approach will depend on the particular problem under consideration. 

The procedures sketched above may be continued to any desired order of approx- 
imation. The approach based on the summation formulas, as in (2.198) and (2.202), 
is a version, in this context, of the standard Rayleigh-Schrddinger perturbation theory. 
The second approach, using (2.204) and solving a differential equation is analogous to 
the Dalgarno-Lewis perturbation theory [60,61] (see also [62,63]). 

However as already mentioned, there are a number of differences between the prob- 
lem of calculating the spectrum of a one-particle operator and that of obtaining the 
spectrum of excited states for a many-particle system. The most obvious is that the 
effective wave equation (2.182) contains the mass operator which is not given a priori 
but is itself defined recursively as the successive terms in (2.183) are obtained. Each 
is expressed in terms of the propagators, or the reduced density matrices, 
p^ . For example, contains the density matrix pf* of (2.186), while contains 
the density matrix 

pf (x u x 2 ,x 3 ,x 4 ) = T \imG { 2 0) (1234), (2.205) 
where the limit is taken such that 

U > h > ti > t 2 . 

For the one-particle propagator, with Fourier transform (2.98), we find 
Go (12) = -tJ2i® (*i2) [l±ra (En)} i 6 (-t 12 ) n(E n } x 

n 

x i/j n (xi) < (x 2 ) exp (-iE n t 12 ) ■ (2.206) 

Substituting (2.206) into (2.158) we obtain the two-particle propagator G 2 \ In view 
of (2.205), this defines the 2-matrix 

p2 ] {xi,x 2 , x 3 , x A ) = n{E m )n (E n ) ip m (x ± ) ip n (x 2 ) x 

m,n 

x [rjxzWM i r m Mr n M] ■ (2.207) 
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Defining the two-particle wave function 

i>mn (xi,X 2 ) := ^= [l/> m (»l) 4>n (x 2 ) ± fa) l\) m (x 2 )} , 

and the pairon density of states 

B mn := n (E m ) n (E n ) , 
the zero-order 2-matrix takes the form 

P2 ] X 2 , X 3 , X 4 ) = B mni>mn {x U X 2 ) ljj* m (x 3 , X A ) . (2.208) 

m,n 

Another difference between the two eigenvalue problems arises from the fact the 
Hamiltonian (2.180) is not necessarily hermitian. Not only can this give rise to com- 
plex eigenvalues signalling the presence of decay, but it also precludes the use of the 
virial or hypervirial theorems which are widely employed in quantum mechanics. These 
theorems may be used in situations, discussed in Section 2.4, in which decay is negli- 
gible and the effective Hamiltonian can be approximated sufficiently accurately by a 
hermitian operator H = K + U . 

If K = —p 2 /2m and p = iV, then if H — K + U is hermitian with eigenfucntions 
■ipn, then the virial theorem asserts that 

The proof of this is based on 

[r • p,H] = 2iK - ir-VU, 

and on the equality 

{ip n , [r • p,H] ip n ) = 0. 

The hypervirial theorem is a generalization of the virial theorem: If if is a hermitian 
operator on Hilbert space with eigenfunctions, ip n , and if A is a well-defined linear 
operator on the same space, then (ip n , [A, H]ip n ) = 0. The proof uses the fact that H 
is hermitian. 

These theorems can sometimes be used to greatly simplify calculations in pertur- 
bation theory for the propagators but only if the effective Hamiltonian is strictly, or 
approximately, hermitian. 

Perhaps the most important peculiarity of the many-body problem is that the 
eigenvalues of the effective equation (2.182) do not give us the spectrum of the excited 
states directly. This spectrum is obtained from the poles of the propagator (2.108) by 
solving the equation 

u n = H n (co n ) - /i (2.209) 

for u n . This defines the single-particle spectrum. Similarly, one can find the poles 
of the two-particle propagator which define the spectrum of the collective excitations. 
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Expanding the two-particle propagator, in analogy with (2.98), one can derive a two- 
particle effective wave equation. One may then proceed to obtain n-particle effective 
wave equations of successive orders n = 1,2,3, .. ., and attempt to solve them by a 
procedure similar to that described for (2.198). Note that the spectra for the effective 
Hamiltonians may depend on thermodynamic variables such as temperature. 



2.9 Screened Potential 

We remarked in Section 2.7, that the unperturbed Hamiltonian (2.152) may include 
potential fields. It is a widely used trick to include in such potentials an averaged form 
of the interactions among the particles. In this way, the convergence properties of the 
perturbation procedure can often be significantly improved. 

An obvious candidate for a potential to be included in the unperturbed Hamiltonian 
for this purpose is the Hartree potential (2.78) which has the form 

V H (1) = / $ (12) p (2) d (2) = ±i J $ (12) G (22) d (2) , (2.210) 

where $(12), introduced in (2.63), is often referred to as the bare interaction potential. 
The unperturbed inverse propagator can then be written as 



Go 1 (12) 



5(12), (2.211) 



in which we have introduced the effective potential 

V(l) :=V H (l)-fi(l). (2.212) 

Thus, already in the zero approximation we partially take account of interaction among 
the particles by means of the self-consistent field (2.210). It is therefore appropriate 
to call this the Hartree or the mean-field approximation. 
The total inverse propagator then takes the form 

G _1 (12) = G^ 1 (12)-S'(12), (2.213) 

in which the supplemental term, 

£' (12) := E (12) - V H (1) 6 (12) , (2.214) 

is called the exchange-correlation self-energy since it takes into account exchange and 
correlation effects not included in the Hartree approximation. We may then proceed 
further by following the perturbation scheme outlined in Section 2.7. 

However, there is an alternative procedure which makes use of the variational deriva- 
tives discussed in Section 2.6. In following this route an important role is played by 
the screened potential 

W(12):= (13)$(32)d(3) , (2.215) 
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in which 

is the so-called inverse dielectric function. The screened potential (2.215) is sometimes 
also called the shielded potential; and the function e _1 (12) in (2.216), the dielectric 
response function. Further, the solution, £(12), when it exists, of the integral equation 



J ^ 1 (13) £ (32)rf(3) =5(12), 



is called the dielectric function. This is an appropriate name since it is a generalization 
of the familiar dielectric constant of electro-magnetic theory. 
To obtain the inverse dielectric function we evaluate 

5V(1) . , io , , . /■_,.,_. 5G(33) 



<5/i(2) v ' J y J 5/2(2) 

and 

5G(12) r 5G (12) (4) 



r OG (12) 01/(4) 
7 5\/(4) 5u(3) 1 ' 



5/x(3) y <SV(4) 5/x(3) 
Introducing the polarization function 

n(i2)-Mi!_ ± ,^(ii) (2217) 
iil ^ j -5y(2)" ± ^y(2)' 

we finally have an integral equation 

e' 1 (12) = 5 (12) + J $ (13) H~ (34) e' 1 (42) d (34) (2.218) 

satisfied by the inverse dielectric function. 

This equation must be supplemented by an equation for the polarization function 
(2.217). Varying the equation of motion (2.69) with respect to the effective potential 
(2.212), we are led to 

/^M G( 3 2)d( 3 )+ / G - 1(1 3)§f d( 3,.0, 
which implies that 

«?(12) t^„,JG- 1 (45). 



6V(3) J v ' 5V (3) 
In order to calculate the polarization function we shall need the vertex function, 

r(123) = -%lT' (2 - 219) 
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which derives its name from the graphical representation of Green's functions by Feyn- 
man diagrams. We find 

n (12) = ±i J G (13) G (41) T (342) d (34) . (2.220) 

In order to make use of this formula we need an equation for the vertex function. A 
variation of (2.213) leads to 

^(12) _ 5m)sm) *S'(12) 
8V(3) ~ 6 WOW S V(3) ' 

and 

SZ'(12) f SSf (12) SG (45) 



5V(3) J <5G(45) 5V(3) 
Using these in (2.219) gives 



Since, as was shown above, 

ft? (12) 
5V{3) 

we finally obtain the integral equation 



J G(14)r(453)G(52)rf(45), 



T (123) = 5 (12) 5 (13) + / 5 ^1 12 } g (46) G (75) V (673) d (4567) (2.221) 

J oG (45) 

for the vertex function. 

We still need an equation for the exchange-correlation self-energy (2.214). Using 
the expression (2.70) for the self-energy and the variational representation (2.150) for 
the two-particle propagator, we have 

£' (12) = -i ( $ (13) 5Gi \^ G~ l (42) d (34) . (2.222) 
J djj, (3) 



By varying (2.69) with respect to /i(l), we find 

/ff G -(3 2 ) d ,3) + / G(1 3)^f) d( 3).0. 

Therefore (2.222) can be rewritten as 

T! (12) = i [ $ (13) G (14) — e (34) . (2.223) 

J dji (3) 
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Using the relation 

SG- 1 (12) _ r 5G- 1 (12) 5V (A) 
5/x(3) ' J 8V{A) 5/2(3) 1 h 
and the definitions (2.216) and (2.219), the equation (2.223) takes the form 

£' (12) = i J G (13) T (324) W (41) d (34) . (2.224) 

We have thus obtained a closed system of five equations (2.215), (2.218), (2.220), 
(2.221), and (2.224) for W,e _1 ,n,r, and respectively. The equations (2.215) and 
(2.218) may be combined to yield the equation 

W (12) = $ (12) + J $ (13) n (34) W (42) d (34) (2.225) 

for the screened potential. 

To provide a conspectus of this set of equations we write them in symboolic form 
as follows: 

W = e- 1 $, 
e" 1 = l + Me- 1 , 
IT = ±iGGT, (2.226) 

r = i + 5 -^ggt, 

S' = iGTW. 

This system of equations may be solved by an iteration procedure, starting from 
the Hartree approximation corresponding to 

So (12) = 0. (2.227) 

Substituting this into the right-hand side of (2.226) gives 

T (123) = 5(12)5(13) (2.228) 

and hence 

n (12) = ±iG(12)G(21). (2.229) 

Let e _1 (12) = 5(12), then 

W (12) = $ (12) . (2.230) 

This implies that screening is absent in the zero approximation and the second equation 
of (2.226) gives the inverse dielectric function 

e~ H l F (12) = 5 (12) + J $ (13) n (32) d (3) (2.231) 
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in the Hartree-Fock approximation in which the exchange-correlation self-energy is 

E' HF (12) = iG (12) $(12). (2.232) 

Since we began with Hartree approximation, the propagators in (2.231) and (2.232) 
are also to be taken in this approximation. 

Substituting (2.229) into the second equation (2.226) gives the equation 

a],], (12) = 5 (12) + J $ (13) II (34) £r % (42) d (34) , (2.233) 

the solution of which is the inverse dielectric function in the random-phase approxima- 
tion. The corresponding screened potential is 

W RP (12) = J e RP (13) $ (32) d (3) . (2.234) 

Inserting this into the fifth equation (2.226) defines the exchange-correlation self- 
energy 

S' 5P (12) = iG{12)W RP (12) (2.235) 

in the screened-potential approximation. 

The iteration procedure can be continued employing either (2.232) or (2.235). How- 
ever, successive approximations soon become overly complicated. In order to avoid 
excessive complications, it is sometimes more convenient to try to include in the zero 
approxiamtion more information about the properites of the system. For example, 
instead of starting with the Hartree approximation (2.227), we might begin with the 
local-density approximation 

Z' LD (12) = V LD (l)5(l2), (2.236) 

in which the local-density potential 

V LD (l)=V LD \p(l)\ (2.237) 

is a functional of the local density p(l) such that (2.237) accounts for part of the ex- 
change and correlation effects. Because of this, (2.237) is also called the local exchange- 
correlation potential. 

Variation of (2.236) with respect to the propagator G gives 

^|1 = ±^(12)K LD (23)5(34), 

where 

Kld (12) := (2.238) 
It then follows by means of (2.226) that we can express the vertex function in the form 

r(123) =5(12)r £D (23), (2.239) 
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where the new function, Tld(12), called the local-density vertex, satisfies the equation 

T LD (12) = 5 (12) + J K LD (13) n (34) T LD (42) d (34) . (2.240) 

Substituting for T(123) in the third equation (2.226) leads to the following formula for 
the polarization function 

U LD (12) = J n (13) T LD (32) d (3) (2.241) 

in the local-density approximation. In this approximation, the inverse dielectric func- 
tion is given by a solution of the integral equation 

e L l D (12) = 5 (12) + J $ (13) n (34) r LD (45)e L ^(52)rf(345). (2.242) 

Consequently, the screened potential in the local-density approximation is 

W LD (12) = J e L l (13) $ (32) d (3) . (2.243) 

Substituting this in (2.226), we define the screened potential 

W LV (12) := J T LD (13)W LD (32)d(3) (2.244) 

in the local-vertex approximation. It follows from (2.240) that the potential (2.244) 
satisfies the equation 

W LV (12) = J Y LD (13) $ (32) d (3) + 

+ J T LD (13) $ (34) ITo (45) W LV (52) d(345). 
With a solution of this equation we may define the local-vertex approximation, 

S' Ly (12) = iG (12) W LV (21), (2.245) 

for the exchange-correlation self-energy [65-69]. 

Since we started from the local-density approximation (2.236), we should take the 
propagators in (2.244) and (2.245) in the same approximation. 

Analysing the structure of perturbation series in quantum electrodynamics, Feyn- 
man suggested a graphycal interpretation of perturbative terms in the form of diagrams. 
Such an interpretation permitted to understand better the physical meaning of pertur- 
bative expansions. It helped also to separate an infinite class of diagrams which could 
be summed, thus, renormalizing and improving perturbation theory. This resummation 
is usually done by selecting those diagrams which correspond to a set of terms forming 
a geometrical progression. Together with the penetration of the field-theory methods 
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into the many-body problem, the diagram technique has become popular among the 
statistical-physics community [70]. 

We shall not consider here such diagrams because of the following: A diagram is 
nothing but the illustration of an analytical expression. It is possible to accomplish 
calculations without diagrams, while, dealing with diagrams, one always has to return 
finally to analytical formulas. In many cases the diagrammatic illustration helps to 
clarify the resummation procedure. However, the latter can always be done without 
invoking the former. 

2.10 Heterophase States 

Till now, speaking about a statistical system we have tacitly assumed that the same 
type of order occurs in the whole volume of the system. It may happen, however, 
that different types of order appear in different parts of the system. If these regions of 
different order are macroscopic then one says that the system is divided into domains. 
In this case, each domain plays the role of a separate system, so that the problem is 
actually reduced to that of a system with one fixed order, with a slight complication 
coming from the boundary conditions describing interdomain layers. 

A much more complex situation develops if the regions of different order are of 
mesoscopic sizes and, moreover, are randomly distributed in space. Since each type of 
order is related to a thermodynamic phase, we may say then that the system is het- 
erophase. The mesoscopic regions of different phases not only are randomly distributed 
in space but may also oscillate in time. Because of this, they are called heterophase 
fluctuations. To describe the system with such fluctuations is a problem far from triv- 
ial. In this section we present a general approach for treating statistical systems with 
heterophase fluctuations [19,71-80]. 

Consider several thermodynamic phases which we enumerate by the index v = 
1,2,.. .. Recall that different thermodynamic phases correspond to qualitatively differ- 
ent statistical states. At the same time, to each phase one can put into correspondence 
a space of typical states which is a Hilbert space Ti u of microscopic states h„ G 7i u such 
that one can define an order operator fj whose mathematical expectation 

(h v ,r)h v ) = r){h v ) 

gives an order parameter r)(h u ) which is qualitatively different from rj(h^) if fj, ^ v. 
The direct sum 

X = ® V U V (2.246) 

describes a system which can be in either one or another thermodynamic phase. In 
order to describe the states of coexisting phases we need to deal with the space 

y = ® v Hv (2.247) 

which is the tenzor product of H u . The spaces X and y are uniquely isomorphic but 
correspond to different physical situations. 
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Let a statistical system in the real-space region V C R 3 be divided into subregions 
Vj, C V filled by different thermodynamic phases. A family {Vj,} of subregions 
forms a covering of the region V when 

V = U„V„, V = Y,V V , (2.248) 

V 

where 

V = mesV, V v = mes\ ' v . 

The space separation of subregions \ u can be characterized by the manifold indi- 
cator functions 

^ : ={i:3v; (2 - 249) 

A set 

£:= {&r)\ re V,i/ = 1,2,...} (2.250) 

of the indicator functions ( |2.249|) uniquely defines a phase configuration in real space. 
Operators of the density observables acting on H v have the general form 

CO , 

x A uk (rj ll r 2 , . . . , r k )d r 1 dr 2 . . .dr k , (2.251) 

in which A vk {. . .) is an operator distribution. Here and in what follows the integrals 
in real space are assumed to be over the whole region V, if not specified. Operators of 
local observables, for a subregion A, are given by 

A V (£,A)= ( A u (tr)dr. (2.252) 

J A 

The operator density Q2.251 ) can be considered as a limiting case of Q2.252 ) in the sense 
of the limit 

mesA^o mesA 

The local operators acting on the space ( |2.247| ) are defined as 

A(£, 7) = ® V A V (£, r ), A(£, A) = @ V A V (£, A). (2.253) 

The family 

A(0 := {A^,r), A(£,A)\ re V, A c V} (2.254) 

of all operators of local observables forms the algebra of local observables. 

When a phase configuration characterized by the set ( |2.250 ) is fixed, then we have 
the case of frozen phase separation. However, when the phase configuration is randomly 
distributed in space and, in addition, may fluctuate in time, we have to consider the set 
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(|2.250|) as a stochastic variable. Then observable quantities are given by the averages 
of operators from the algebra of local observables ( |2.254|) , with the averaging procedure 
containing the trace over microscopic degrees of freedom and, also, an averaging over 
phase configurations. The latter is a functional integral over the manifold indicator 
functions Q2.249| ), which we shall denote by / Z>£. The statistical averaging, as usual, 
must include a statistical operator p(£). Thus, the observable quantities are given by 
the average 

(A) =TrJ p(£)A(Z)V£. (2.255) 

Here and everywhere in what follows it is assumed that the operators of taking trace 
and of integrating over £ always commute. 

As examples of the operators ( \Z.'2b'd ), we may write the operator of the energy 
density 

£(£,r) r), (2.256) 

where 

E„(£, 7) = it (£, 7) + V„(Z, 7) (2.257) 

is the operator of the energy density for a i/-phase, consisting of the kinetic-energy 
operator density. 



and of the potential-energy operator density 



2 



Mr) (2.258) 



m,7) = l -J Ur)Ur'Wu{rW v {7')x 

x $(? - 7')ip u (7')ip u (7)d 7' . (2.259) 
Another example is the operator of particle density 

N(Z,7)=® U N U (Z,7) (2.260) 

with 

AUC, r) = Ur)4(r)Mr )■ (2.261) 

To calculate the average ( |2.255| ), we need to know the statistical operator p(£). 
This can be defined by invoking the principle of unbiased guess, that is, by maximizing 
the information entropy. The entropy itself is defined by the expression 

S(p) = -Trj p(£)lnp(£)P£. (2.262) 

As additional conditions for p(£), we require that this statistical operator is normalized 
so that 

Tr J p{i)Vi = 1, (2.263) 
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that the local energy density is given by 



E(r)=Tr J P (£)E&t)V£, 



and that the local particle density is 

N(r) = Tr \ p{£)N{^ r)D£. 
For the informational entropy we have 



(2.264) 



(2.265) 



S inf = S(p) + c[TrJp(OV^-l 
+ J /3(r) \E(r) -Tr J r)V£ 



+ 



d r + 



+ / l{r) 



N(?)-Tr p(ON(Z,r)n 



d r. 



(2.266) 



where (, (3(r), and 7(r) are the Lagrange multipliers. Conditions ( |2.263| ), (|2.264j) , 
and ( p. 265] ) are recovered from the variational equations 



dS 



inf 



5 Si 



inf 



5S 



inf 



0. 



The statistical operator p(£) is defined as a maximizer of ( 2.266|) , that is, from the 
equations 

oSinf „ S Sj /r , 



'inf 

MO 



o. 



'inf 

¥(0 



< o. 



(2.267) 



Introducing the notation 

C = l-lnZ, 7 (f) = -(3{r)n{? 
from the extremum condition of ( |2.267| ) we obtain 



where 



is the partition function and 



P(0 = ^exp{-Q(0}, 



Z = TrJ exp{-Q(£)}D£ 



Q(0= P(r) E(£,r) - p(r)N{£,?) dr 



(2.268) 
(2.269) 

(2.270) 



is the quasi-Hamiltonian. 
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Now we need to concretize the procedure of averaging over phase configurations, 
that is, we have to define explicitly the functional integration over the stochastic vari- 
able ( |2.250| ). All conceivable variants of the phase configurations form a topological 
space /C = {£}. To integrate over the stochastic variable £, we have to define a func- 
tional measure / on the space /C. 

Note that averaging over space configurations contains two types of actions. The 
first one averages over all possible phase configurations under a fixed set 



of geometric probabilities 
with the evident properties 



x := {x u \v = 1,2,...} (2.271) 



Xv := V I u7)d 7= v (2 - 272) 



5^ = 1, 0<x v <l. (2.273) 



The second action is the variation of each probability ( [2.272 ) from zero to unity taking 
account of their normalization ( |2.273| ). 

The resulting differential functional measure T>^ can be written as 

Vi = dxV x i (2.274) 

with 

dx = 6 (Y, x » ~ lj II dx »i (2.275) 

where x u G [0, 1]. 

To define T> x ^, let us divide the whole system into n parts, so that each subregion 
V„ is covered by a family of n^-subregions \ ui such that 

n v n v 
i=l i=l 

and with 

n = ^n u , v ui = mesV vi . 

V 

The shape of each sub-region V ui can be arbitrary. For each V„j, let us fix a vector 
a^iG V ui defined as the center of the local coordinate system in V^. Then the manifold 



indicator function ( |2.249| ) can be presented as 

Ur) = Y,Z*(r-a vi ), (2.276) 

i=l 

where 

'1, re V„i 

£ ui (r - a ui ) = < 

, 0, rgV*. 
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It is clear that if the sub-subregions are small enough, then they can be used as 
blocks for producing any phase configuration. A displacement of the blocks can be 
accomplished by moving the block centers a ui . 

In this way, we come to the definition of the functional differential 

^:=nftv (n-°°), (2-277) 

V i=l ' 

in which fl„6 V and the limit n —>■ oo implies that 

n — ► oo : n u — > oo, i^i — > 0, — > const. (2.278) 



The limit (|2.278|) is to be taken after the corresponding functional integration. Thus, 



/ mv*t = ii- / m n n (2-279) 



for a function of £ denoted by we have 

n 

^ i=i 

For example, it is straightforward to check that 

/ V x £ =1, J U rYPX 

More generally, the following statement holds. 
Proposition 2.1. 

Consider a class of functionals of the form 

oo - 

F (0 = £££•••£ U (ri)^(r 2 ) . . . ^(r fe ) 

fc = fl ^2 ffc 



X 



x F UlU2 .. Mk (r 1 ,r 2 ,...,7 k )dr 1 dr 2 ...dr k . (2.280) 
Then the functional integration ([2.279 ) gives 



jF(£)V a> Z = F(x), (2.281) 



where F(x) has the form of fl2.280| ) but with all £v{ r ) changed by x v . 

Proof is based on a direct integration of products of the manifold indicator func- 
tions (ET276D over the block vectors a ui with the differential measure (|2.277|) . In this 



integration, the products of ^ with coinciding vectors a vi yield Lesbegue zeros. For 
example, 



£ J Z,vi{r\ - a ui )£ ui (r 2 - a ui )V x £ 
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lim ^ / ^„(ri - a^)^(r 2 - a^)— = x„5 nr2 , 



where 5 rir2 is the Kronecker delta. Therefore, only the manifold indicator functions 
with different block vectors a ui contribute to ( |2.279| ) which leads to ( |2.281| ). 



Proposition 2.2. 



The partition function ( |2.269| ) with the quasi-Hamiltonian ( |2.270| ) and the differ- 



ential functional measure ( p.274 ) can be written as 



Z = Tr J exp{-Q(x)}dx, (2.282) 



where Q(x) has the form of ( |2.270| ) but with all £ v (r) changed by x V} and dx is defined 
in ( j2~275l ). 



Proof involves the previous proposition after expanding the exponential in (|2.269| ) 
in powers of Q(£) and obtaining 



J exp{-Q(£)}XU = exp{-Q(x)}. 



The integration in ( |2.282| ) with the differential measure ( |2.275| ), if accomplished in 
the thermodynamic limit iV — ■> oo, defines a natural thermodynamic potential for a 
heterophase system as follows. 



Proposition 2.3. 

Let a function 



y(x):=- lim — InTr exp{-Q(x)} (2.283) 



of x in the hypercube 

H := {x\ J2 X » = !> < ^ < 1} 



have an absolute minimum 



y{w) = absminy(x), (2.284) 



in which the set 

w:={w v \u = l,2...} (2.285) 
consists of elements with the property 

^2w u = l, 0<w v <l. (2.286) 

V 

Then, in the thermodynamic limit, 

y(w) = - lim i-lnZ. (2.287) 

JV — >oo _/V 
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Proof uses the asymptotic, as N — > oo, presentation 



exp{— Ny(x)}dx 



for the partition function ( J2.282|) with definition (2.283). Then, the Laplace method 
of integration immediately yields (|2.287|) . 

Thus, we come to the conclusion that the role of a natural thermodynamic potential 
for a heterophase system is played by the function 



y(w) — — lim — InTr exp{— Q(w)}. 



(2.288) 



The set (|2.285|) consists of the geometric probabilities w v defining relative weights of 
each thermodynamic phase. Consequently, w u may be called the phase probability. To 
find these probabilities, we have to minimize the heterophase thermodynamic potential 
( |2.288D . The minimization means that the probabilities w v are to be defined by the 
equations 

d 2 
-y(w) = 0, — 2/M > (2.289) 



dw v K 



under condition ( |2.286|) . In addition, one has to analyse the thermodynamic potential 
y(w) for w on the boundary of the hupercube H, that is, when one or several w u are 
equal to one or zero. The first equation in (|2.289| ) can be called the equation for phase 
probabilities, and the inequality in ( 2.289|) is the condition of heterophase stability. 

Taking account of the form of the heterophase thermodynamic potential (|2.288| ), 
as the equation for phase probabilities we get 



lim 



1 







■QH 



N^oo N \ dw u 
And the condition of heterophase stability becomes 



(2.290) 



lim — 

N^oo N 



' d 2 



:Q(w) 



' d 
dw v 



Q(w) 



> 0. 



(2.291) 



The quasi-Hamiltonian Q(w) is defined in ( [2.270| ) which together with ( p.256|) -( p.261 
gives 

Q(w) = ® v Q v {w), (2.292) 



where 



2m 



+ U(r)- fi(r) 



ip u (r)d 7 + 



+ ]zwl J (3(r )i)t(r )ipt(r )$(? - r')ip u (r')ilj lJ (r)d r dr' 



(2.293) 



Finally, we need to find an effective presentation for the expectation values ( |2.255| ) 
of the algebra of local observables ( p.254| ). 
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Proposition 2.4. 

Let the conditions of the Proposition 2.3 be valid, then the expectation values 
Q2.255Q asymptotically, as iV — > oo, take the form 



(A) Trp eff (w)A{w) (N -> oo), (2.294) 



where A(w) is the algebra of local observables Q2.254D with all £ v (r) changed by w U) 
and the effective statistical operator is 

PeffM := J—exp{-Q(w)} (2.295) 

with the effective partition function 

Z eff :=Trexp{-Q{w)}. (2.296) 

Proof starts from Proposition 2.1, according to which 

' p(t)A(0T> x { = P (x)A(x), 



where ^ 

pO) = ^exp{-Q(a;)} 



with Z from ([2.282 ). Then, introducing 



._ Trexp{-Q(x)M(x) 

exp{-A^)} ' {Z - Z3<) 



we have ^ 

Trp(x)A(x) = — exp{— Ny(x)}A(x). 

Zj 

Using the method of steepest descent, for iV — » oo, we obtain 
exp{— Ny(x)}A(x)dx ~ exp{— Ny(w)}A(w) 



27T \ 1/2 



where the product over v contains the number of factors by one less than the number of 
phases, since one of w v must be expressed through others because of the normalization 
( gg86[) , and 

„ _ d 2 y{w) 
dw 2 



y" - -Br > °- 



Similarly, we find 

Z ~ exp{-Ny(w)}]J 
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It follows that 



1 

~Z 



exp{— Ny(x)}A(x)dx ~ A(w) 



as N — > oo, that is 

(A) ~ A(w). 

Remembering the definition Q2j97| ), we come to ( gggj ) with (|2~295l) and ( g29g ). 

When there are no external fields specially inducing nonuniformity, the local inverse 
temperature and chemical potential are to be uniform through the system: 



(2.298) 



These equations have the standard form of equilibrium conditions. However, we should 
not forget that a heterophase system, as is explained above, is quasiequilibrium. Condi- 
tion (|2.298|) appears after averaging over phase configurations. Therefore this condition 



has the meaning of the condition of equilibrium on average or it may be called the con- 
dition of heterophase equilibrium. 

Under condition ( |2.298|) , the quasi-Hamiltonian ( [2.292|) becomes 



Q(w) = (3H{w) 
with the renormalized or effective heterophase Hamiltonian 

H{w) := @vH v {w) 

consisting of the terms 



(2.299) 



(2.300) 



HJw) 



V 2 

-— + U(r 
2m 



i) u (r)d r + 



+ ^wl J i[)l(r)i{;t(r')^(r - ~r')i) u {r')^ u (r)d r dr' 



(2.301) 



Emphasize that ( p.300| ) represents not just one system with a given phase separation 
but an infinite number of systems with all possible phase configurations. Because of 
this, each term ( 2.301|) may be called the phase-replica Hamiltonian. Owing to the non- 
linear renormalization in ( |2.301| ), interphase effects are included into the description. 
The interphase state is defined as 



(A) 



hit 



(A) - J2w u (A) u , 



(2.302) 



that is, as the difference between the heterophase state (^4), and the respectively 
weighted sum of pure states 



(A)u 



lim (A). 
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As an example, consider two coexisting phases. Then the phase-probability equa- 
tion (|2.290|) , with the quasi-hamiltonian (|2.299|) gives 



Wi 



2$ 2 + K 2 -K 1 + v{Ri ~ Ba) 



2($x + $ 2 



w 2 



1 - Wu 



(2.303) 



where 



K " = N 



V 2 

2m v ' 



V>„(r) )d r. 



Thus, all phase probabilities are defined in a self-consistent way. 

Each field operator ip u (r) is defined on a weighted Hilbert space Ti v and has the 
standard commutation or anticommutation relations, 



0. 



depending on the kind of statistics. Field operators defined on different spaces, by 
definition, commute with each other: 







^/v(V),Vv(r ) 
Therefore, the evolution equation for 

Mr,t) = U+m v (r)U v {t) 
has the usual Heisenberg form 





iy ± i/) 



(2.304) 



U v )t) 



-iH v t 



For each phase replica, indexed by v = 1, 2, . . ., we may introduce a propagator 

G„(r, t- r ', f ) := -i(f^„(r, t)^ v (r , t')) . (2.305) 

Similarly, we may introduce many-particle propagators. Then we can write the evolu- 
tion equation for ( |2.305|) and proceed in analogy with one of the techniques described 
in the previous sections. Dealing with the effective Hamiltonian ( |2.30C| ), we have the 
dynamical equations for different phase replicas formally decoupled, although thermo- 
dynamically all equations are related with each other through the geometric proba- 
bilities (|2.285|). The possibility of realizing this dynamical decoupling is one of the 



main advantages of the approach formulated in the theory of heterophase fluctuations 
[19,71-80]. 
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3 Correlated Iteration 



The Green function technique expounded in the previous Chapter can be applied to 
many statistical systems. However, there exists an important class of systems which 
exhibits such a strong correlation between particles that a direct application of the 
methods discussed in the previous Chapter becomes impossible. For such systems 
there should be developed a special approach taking into account strong interparticle 
correlations from the very beginning. 

In this Chapter, we shall present such an approach to treating strongly correlated 
systems, based on Refs. [81-84]. 

3.1 Singular Potentials 

The iteration procedures for developing a perturbation series which we have discussed 
hitherto are convenient for systems of particles interacting through non-singular or 
weakly singular potentials. We make these terms more precise by classifying potentials 
according to their behaviour as r\ 2 — * 0, where r\ 2 = | r 12 1 , and ri 2 = r 2 — X\. We shall 
say that $>(xx,x 2 ) is a non-singular potential if 

lim |$(xi,x 2 )| < 00. (3.1) 
A singular potential is one for which 

lim |$(xi,x 2 )| = 00. (3.2) 
>~i2^o 

A weakly singular potential is singular but is integrable in the sense 



n 



< 00 , (3.3) 



where TZ is a finite region including r 12 = 0. A strongly singular potential is singular 
and is not integrable so that 



$(x u x 2 )dr l2 



00. (3.4) 



For example, if the behaviour of a potential at short distances is governed by the 
power law 

$(xi,x 2 )~— (ri2->0), 

and TZ is a small finite ball centered at the origin, then the potential, $ is nonsingular 
if n = 0, weakly singular if < n < 2, and strongly singular if n > 2. Thus the 
Coulomb-like potentials with n = 1, are weakly singular. Whereas the Lennard- Jones 
potential with n — 12, which is sometimes used to model interaction of molecules, is 
strongly singular. 
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For systems interacting via a strongly singular potential, any iteration procedure 
starting from the Hartree-Fock approximation breaks down since in the very first step 
the Hartree potential is divergent. This is because the input to the Hartree approx- 
imation assumes independent particles and does not take into account short-range 
correlations which could effectively shield the singularity. These short-range correla- 
tions might be retained by a decoupling procedure in which the three-particle propa- 
gator is approximated by a linear combination of two-particle propagators. However, 
decoupling at any stage effectively precludes subsequent approximations. 

What we really need is a perturbation procedure which satisfies the following con- 
ditions: 1) it should start from an approximation containing no divergencies and must 
therefore take into account strong short-range correlations, 2) it should provide an algo- 
rithm which permits recursive calculation of higher-order corrections. We now propose 
a procedure which we believe meets these desiderata and which we have named the 
correlated iteration theory. 

3.2 Iterative Algorithm 

Let us consider a situation with a strongly singular interaction potential $(12) and let 
be given by (2.75). Suppose we can find a function, s(12), such that 

GT(1221) := -s(12)p(l)p(2), h > t 2 , (3.5) 

is an approximation for the 2-particle propagator. We shall call s(12) a smoothing 
function. It follows from this equation that 

s(12) = s(21). (3.6) 

We shall require that the smoothed potential 

$(12) := s(12)$(12) (3.7) 

be nonsingular or weakly singular and therefore integrable in the sense of ( p.3|) . 
We now introduce a function, .0(123), satisfying the relation 

G 2 (1223) := s(12) / L>(124)G(43)d(4), (3.8) 

which we call the doubling function. This function will be used to transform a 1- 
particle into a 2-particle propagator, while the smoothing function takes inter-particle 
correlations into account. An initial approximation, D , can be obtained by finding a 
solution of the integral equation obtained by combining ( |3.5|) with ( |3.8| ), giving 

G°(1223) = J D (124)G(43)d(4). (3.9) 
It is not difficult to verify that 

D (123) := 5(13)G(22) ± G(12)5(23) (3.10) 
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is a solution of (|3.9|). Using ( ft.8|) , the exact self-energy becomes 



£(12) = ±i /" $(13)D(132)d(3) (3.11) 

which, in view of ( |3.10|) , has as first approximation 

E x (12) = z$(12)G(12) + 5(12) J $(13)p(3)d(3). (3.12) 

This differs from the Hartree-Fock form in that the bare interaction potential $ has 
been replaced by the smoothed potential Therefore, no divergencies appear in ( |3.12| ) 
which may be called the first correlated approximation or the correlated Hartree-Fock 
approximation. 

In order to define a recursive procedure for obtaining higher-order approximations 
we shall need a set of equations which determine the doubling function and the self- 
energy. Establishing these equations is the next, and most complicated step in our 
theory. By a variation of the equation of motion (2.69) with respect to /i, we have 

SG(12) r 5G- 1 (45) 



5[i(3) J 5/i(3) 
From the definition of the inverse propagator it follows that 

SG-\12) , n9Um , g(12) 

where 

^(12) _ r jglg) 8G(45) 
5(i(3) J 5G(A5) 5fi(3) 1 '' 

Therefore, 

<?G(12) 



M3) G < 13 > G < 32 > + 

Substituting in this the expression 

^^ = ±[G(12)G(33)-G 2 (1332)], 

which results from the variational representation (2.151), we obtain the equation 

G 2 (1223) = G£(1223) + (3.13) 

+ / G ( U )^^ [C 2 (6227) - G(67)G(22)] G(53)rf(4567) 
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for the 2-propagator. By using the relation (|3.8| ) which defines the doubling function, 

we get 

s(12)D(123) = D (123)+ (3.14) 

<5S(43) 
'<JG(56) 

This suggests that we should define an integral operator X on H 3 as follows 



+ / G(14)^|[ s (52)D(527)G(76)-G(22)G(56)5(67)]rf(4567). 



X/(123) := s(12)/(123)+ (3.15) 

+ / G(14)^|^[G(22)G(56)<5(67)- S (52)/(527)G(76)]rf(4567). 
Then ( p.14 ) may be rewritten as 

XD(123) = D (123). (3.16) 
Now define, as left inverse of X, the operator Y such that 

Y:=X _1 , YX = 1, (3.17) 
where 1 is the identity operator. Acting with Y on ( |3.16| ), we obtain 

£>(123) = Y-D (123), (3.18) 

which is the main equation for the doubling function. 

To obtain an iteration procedure for solving ( 3.18| ), we proceed as follows. Assume 
that our initial approximation, (|3.5|) , models the major properties of the system under 
consideration so that D will be close to D . In this case, we shall assume that X is 
close to the identity. We therefore try to obtain Y by expressing it as 

Y=[l-(l-X)]-\ (3.19) 

which leads to the Neumann series expansion 

oo 

Y = ^(l-X) fc . (3.20) 

fc=0 

From (p. 15 ) we see that X is a functional of the self-energy S, 

X = X[£]. (3.21) 

Let 

X fe = X[£ fc ] (X := 1) (3.22) 

be the same functional acting on the fc-approximation to the self-energy. Define the 
^-approximation to Y by 

k 

Y fc :=£(l-X fc )P (Y :=l). (3.23) 

p=0 
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If the sequence {Y fc } has a limit, then 



lim Y fc , (3.24) 

k— >oo 



provided that E = lim^oo 

With equations ( |3.11|) , 03.18Q , and (|3.20|) we have reached our immediate goal of 
being able to define an iteration process for the calculation of the doubling function 
and the self-energy. The recursion can be organized according to the scheme: 

Dk —> — > Yfc+i — > -Dfe+i. (3.25) 



For k = we start from ( |3.10|) leading to ( |3.12| ). Then because of the variational 
derivative 

= ,1,(14) [5(13)5(24) ± 5(12)5(34)] , 

we find 

Y x =2-Xi. (3.26) 



It follows that (|3.18|) gives 

£>i(123) = [2 - s(12)]D (123)± 
±i J [G(14)G(23) ± G(13)G(24)] G(42)<l(43)d(4)T (3.27) 
Ti I [G(13)G(44) ± G(14)G(43)]G(22) [1 - s(42)] $(43)d(4). 



Combining this with (|3.11|) yields the second correlated approximation for the self- 
energy: 

E 2 (12) = z$(12)G(12) + 5(12) / $ 2 (12)p(3)d(3)- 



- J $(13)G°(1324)G(43)$(42)rf(34)+ (3.28) 

+ J $(134)G 2 1 (1332)G'(44)$(32)d(34), 
which involves the following two new functions, the twice smoothed potential 

$ 2 (12) := $(12) [2 - s(12)] = s(12)[2 - s(12)]$(12), 
and the effective triple potential 

$(123) := <5(13)[1 - s(23)] = s(13)[l - s(23)]$(13). 

Provided that the smoothing function, s(12), is known we may find higher-order ap- 
proximations by following the scheme (|3.25|) . 
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3.3 Smoothing Function 

According to ( |3.5| ) the smoothing function is 



s(l2) = (h 2 > 0). (3.29) 

In view of the relations 

G(ll) = =fip(l), G(12)G(21) = ^n (12), 

the denominator in (|3.29|) is given by 

p(l)p(2) = -G°(1221)-ffl (12). 

Since it is clear that the diagonal two-particle propagator, G 2 (1221), plays a special role 
in our correlated approximation procedure, we analyze some of its general properties. 
Consider the density-of-particle operator 

n(l) := ^(1)^(1), (3.30) 

which has expected value 

p(l) = (n(l)>, (3.31) 
equal to the density-of-particles. It is easy to check that when t\ = t 2 , 

[n(l),ra(2)] = 0. 

Thus we may define the chronological product 

fn(l)n(2) := ( } > \ 2 . (3.32) 

w w | n(2)n(l), h <t 2 v 1 

We may then, for the diagonal propagator G 2 , obtain the expression 

G 2 (1221) = -(fn(l)n(2)}. (3.33) 

Note the symmetry property 

G 2 (1221) = G 2 (2112), (3.34) 

which holds both for t x > t 2 and ti < t 2 - This is in contrast to the density-density 
correlation function 

#(12) := (n(l)n(2)) (3.35) 

for which, in general, R(12) = R(21) only when t\ = t 2 . 

The correlation function ( |3.35|) is related to the second-order density matrix by the 
relation 

lim #(12) = 5(xi -x 2 )p(l) + p 2 (x 1 ,x 2 ,x 1 ,x 2 ,t 1 ), 
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which implies that the propagator ( |3.33| ) may be expressed in the form 



G 2 (1221) = -9(t 12 )i?(12) - Q(t 21 )R(21), (3.36) 

where t i2 := t± — t 2 . 

A quantity which is often considered is the density-density Green's function defined 

by 

G p (12) := -i{fn(l)n(2)), (3.37) 

which satisfies spectral representations and dispersion relations similar to those of the 
Bose propagators studied in Sections 2.2 and 2.3. The propagators ( |3.36|) and (|3.37| ) 
differ only by a factor i, 

G p {12) =zG 2 (1221). 
Another important quantity is the pair correlation function, 

J(12) - wmm- ( ' 



Since this is related to ( |3.35|) by the equation 

R(12)=p(l)p(2)g(12), (3.39) 

it follows that g(12) and -R(12) have the same symmetry properties. In particular, 

c/(12) = </(21), (ti = t 2 ). 

By ( |3.38| ), the propagator ( |3.33j ) can be expressed by the formula 

G 2 (1221) = -p(l)p(2) [6(t 12 )0(12) + Q{t 21 )g{2\)} . 

The equations of motion imply that G 2 differs from G 2 only in that it takes into 
account interactions between particles. Thus if $(12) approaches zero, then G 2 (1234) 
approaches G 2 (1234). This implies that the smoothing function ( |3.29j ) must satisfy the 
asymptotic condition 

s(12) -> 1, |xi-x 2 |->oo. (3.40) 

The Bethe-Salpeter equation (2.86) will sometimes enable us to find an approximate 
G 2 PP needed for the equation (|3.29| ) defining the smoothing function. 

3.4 Response Functions 

The functions discussed in this section enable us to characterize the changes in a 
system under the action of some external source or when the corresponding macroscopic 
quantity is varied. For example, the density response function 



119 



describes the change in the average density 

5p(l) = - [ x(12)tf/*(2) 



when the effective chemical potential is varied. 

It is also convenient to introduce the generalized response function 

*< 128 >==*W (3 ' 42) 



of which ( |3.41 ) is the particular case 

X(12)=x(112). (3.43) 
From the variational representation (2.151) we have 

X(123) = % [G 2 (1332) - G(12)G(33)] . (3.44) 
By substituting G\ we find 

Xo(123) = i [G 2 , (1332) - G(12)G(33)] , 
which defines in zero-approximation 

Xo(123) = ±iG(13)G(32) (3.45) 

for the response function (|3.42|) . 

Because of ( |3.43j ) the density response function ( |3.41| ) can be expressed by 



x(12)=i[G 2 (1221) + p(l)p(2)]. (3.46) 
In view of ( |3.34|) this shows that the response function is also symmetric: 



X (12) = X (21) = . (3.47) 

An equation for the determination of the generalized response function can be 
obtained by using ([3.44 ) and ( 3.14 ), giving 



X(123) = xo(123) + J G(14)G(52)^||^| X (673)rf(4567). (3.48) 

In order to make effective use of this equation we need an approximation for S. It 
is perhaps worth stressing again that in the case of strongly singular potentials we 
cannot employ the Hartree-Fock approximation for the self-energy. Some modified 
approximation such as that provided by the correlated iteration procedure is essential. 
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From the correlated Hartree approximation 

S C h(12) := 5(12) J $(13)p(3)d(3), (3.49) 



since 

^ = ^(12)*(23)*(34), 



noting ( |3.48| ) we conclude that 

X(123) = xo(123) + / Xo(124)$(45)x(553)rf(45). 



By specialization this shows that the equation for ( p.43 ) takes the form 



X(12) = Xo(12) + J Xo(13)<l(34)x(42)rf(34), (3.50) 

where 

Xo(12) := ±zG(12)G(21) = n (12). 

We shall say that the solution of ( |3.50|) is the response function in correlated random- 
phase approximation. In the standard random-phase approximation, the smoothed 
potential, would be replaced by the bare potential, If we employ the correlated 
Hartree-Fock approximation ( |3.12| ) for X in ( |3.48| ) we are led to the equation 

X(123) = Xo(123)± (3.51) 
±i f G(14)$(45)[G(42)x(553) ±G(52)x(453)]ci(45). 



We note in Chapter 2 that the poles of the two-particle propagators define the 
spectrum of collective excitations. The relation ( |3.46| ) shows that these propagators and 
the density response functions have the same poles. Since the equations for the response 
functions are usually simpler than those for the 2-propagators, the response functions 
are commonly considered when the collective excitations spectrum is of interest. We 
describe two possible approaches to finding this spectrum. 

1. Substitute an approximation for the self-energy, S, in ( |3.48| ) which results in 
an equation for the generalized response function. If possible, specialize to obtain an 
equation for the density response function. Except for usual cases, this specialization 
is feasible only when the Hartree approximation is employed for the self-energy, but 
then as follows from ( |3.51|) this specialization will be impossible in the next approxima- 
tion. Further, to proceed from the first approximation one needs to solve the resulting 
equation for the density response function. Because of the product in the integrand 
of the right-hand side of ( |3.50| ) this can be quite difficult unless one uses merely a 
random-phase type of approximation. In this latter case we find a solution for x(12) 
in the form 
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where the subscript, CRP, means correlated random-phase approximation. The spec- 
trum of collective excitations is then approximated by the zeros of the equation 



1 - Xo$ = 0. (3.53) 



However, it is not always easy to solve ( |3.51|) for the density response function. 



2. A second approach to finding the spectrum of the collective excitations begins 
by seeking an approximate iterative solution of ( |3.48| ) in the form of a series in powers 
of $. Since the poles of \ are the zeros of \~ x , the spectrum is characterized by the 
equation = after \ is replaced by a polynomial in 5> of degree equal to the order 
of approximation one is using for \. For example, taking 

X - Xo + Xo&Xo (3.54) 

leads to 

X" 1 ^ fl - Xo$) Xo\ (3-55) 



when the infinite series for x _1 is truncated at the first order in <£>. The zeros of 
this equation are the same as those of equation ( |3.53j ). Notice that it might happen 
that ( |3.54| ) has no poles even though ( |3.55|) has zeros! This could result because the 



truncation meant that these two equations are not strictly equivalent. 

In order to carry this second approach through effectively we need to have at our 
disposal a "regular" perturbation procedure, - that is a procedure which does not 
terminate abruptly due to the occurrence of singularities. We can define such a proce- 
dure by employing the same trick for inverting operators as we used in our correlated 
iteration theory for the doubling function and self-energy in Section 3.2. 

We introduce an integral operator, Z, on if 3 such that Zf is defined by 

Z/(123) := / G(14)G(52)^|/(673)rf(4567). (3.56) 



Then the basic equation, ( |3.48 ), for the generalized response function can be written 



as 



(I-Z)x(123) = xo(123). (3.57) 

This is equivalent to 

X(123) = (I-Z)- 1 X0 (123), (3.58) 
in which the operator can be expanded in the series, 



oc 



(I - Z)" 1 = ^ Z*. (3.59) 



fc=0 



It is clear from ( |3.56|) that Z is a functional of the self-energy. We symbolize this by 

Z = Z[S], (3.60) 
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and denote the operator fc-approximation to Z by 



Z fc :=Z[£] fc , Z :=I, (3.61) 

where denotes the ^-approximation for the self-energy in the correlated iterative 
scheme. The ^-approximation for the generalized response function is then Xk, defined 
by 

k 

Xfc (123):=£Z£ Xo (123). (3.62) 

In the previous Section we showed how to obtain successive approximations to the 
self-energy by means of correlated iterations, so following the scheme 

S fc -> Z fc — > Xk (3.63) 

we may obtain successive approximations to the generalized response function. Ac- 
cording to ( |3.43| ) the fc-approximation to the density response function is obtained by 
the specialization 

Xfc (12)=Xfc(H2). (3.64) 
For example, with k = 1 in (|3.62|) we have 



Xi(123)= X o(123)± 
±i I G(14)$(45)[G(42)xo(553) ±G(52)xo(453)]c/(45). 



It follows that 



Xi(12) = Xo(12) + J X o(13)<l(34)xo(42)d(34)± (3.65) 
± / xo(123)$(34) X o(214)rf(34). 



The schemes ( |3.25|) and ( |3.63|) provide a recursive procedure for obtaining a se- 



quence of successive approximations to the density response function in powers of $ 
which define x asymptotically. The poles of the function so-defined, or the zeros of 
characterize the spectrum of collective excitations of the system. 

3.5 Correlation Function 

From the above analysis it follows that the smoothing function and the pair correlation 
function are intimately related with each other. Really, comparing the definitions (3.5) 
and (3.38), we see that the smoothing function can be identified with an approximation 
for the pair correlation function. Therefore, to make the correlated iteration theory 
complete, we need to present a way for finding approximate correlation functions. 
The latter, according to ( |3.38| ), can be defined as 



, , G 2 (1221) . , 

9(12) = -tkW (tl2 > 0) - (3 - 66) 
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If we substitute in (|3.66|) 

G°(1221) = -p(l)p(2)±G(12)G(21), 

then we get 

a(l2)-lT G{12)G{21) (^67) 
flb(12)-l=F p{1)p{2) ■ (3.67) 

A general equation for g(12) can be obtained employing in ( |3.13| ) definition ( |3.66|) , 
which gives 

G(13)G(41) 5S(34) 



#(12) = # (12) 



p(l)p(2) 5G(56) 



x 



G(56)G(22) -#(52) / D(527)G(76)d(7) 



d(3456), (3.68) 



keeping in mind the condition t 12 > 0. Using in fl3.68| ) approximate expressions for £ 
and D, we then come to a set of approximate equations defining #(12). For instance, 
putting S = So = in Q3.68 ) immediately returns us to ( p. 67 ). If we take for £ and D 
the first approximations in ( ft .12] ) and (|3.37|) , respectively, then ( |3.68| ) transforms to 

#(12) = ,/o(12) - J A(1234)#(34)$(34)d(34) , (3.69) 

where we keep only the terms up to first order in <3> and 

A(1234) := i G ^^ {^(32) [2 - #(32)] G°(3224) - G(34)G(22)} ± 

±i G{ ph)^2) ] i 9{42) [2 ~ 9{A2)] G ° (4224) " G ( 44 ) G ( 22 )} • 



Equation (|3.69| ) defines first-order approximation for #(12). 

Let us notice that ( |3.69| ) cannot be iterated starting with #o(12) in the integrand, 
since then integral would diverge because of singularity in $(34). A simplification is 
possible only for A(1234). If we put in the latter # = 1, then 

Taking # = 1 implies that the Hartree approximation is used for which 

G°(1234) G(14)G(23), # (12) 1 . 
With this in mind, we can make a further simplification getting 

1 ^ P(l)p(2) 
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If we substitute in ( |3.68| ) the second-order approximations for £2 and D2, then 



we come to an equation defining the second-order approximation for the correlation 
function g(12). However, this equation is rather cumbersome. Instead, we may improve 
the first-order approximation for g(12) in the following way. We can consider ( 3.69Q 



as an asymptotic expansion in powers of $. An effective limit of an asymptotic series 
can be defined by using the self-similar approximation theory [85-97]. Employing the 
algebraic variant of this theory [98-100], we obtain 

0(12) =So(12)e X p{- J A(1234) -^$(34)^(34) j . (3.70) 

This equation, as is explained above, because of divergences, cannot be iterated 
starting with g Q (12). One direct option would be to solve ( |3.70|) numerically for a given 
interaction potential $. Another way can be by iterating the right-hand side of ( |3.70| ) 
starting with an approximate g app (12) taking into account interparticle correlations. 
For example, we may start [82] with 

g app (12) = \^>(12)\ 2 , 

where \1/(12) is a pair wave function normalized so that to satisfy the condition 

g app {12) -> 1, \x 1 - x 2 \ -> 00 , 

similar to ( |3.40| ). In this way, the correlated iteration theory is completed. 
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Appendix 



It was mentioned in Preface that statistical mechanics is the most general physical 
theory. Here we concretize this statement. 

From the mathematical point of view, one theory is more general than another if 
the latter can be derived from the former by means of some mathematical operations, 
but not conversely. A particular theory can be considered as a limiting case of a more 
general theory. The corresponding limit can be defined explicitly in an unambiguous 
mathematical sense. The procedure of taking a limit is, of course, to be accompanied 
by the related physical interpretation. 

The classification of physical theories can be done by considering three parameters: 
velocity of light, c, Plank constant, h, and Boltzmann constant, k B . The choice of 
these parameters is caused by the fact that each of them is related to a principal 
physical notion. The light velocity, according to the Einstein postulate, is the same in 
all reference systems. The Plank constant defines a minimal quant of action. And the 
Boltzmann constant connects the entropy with a statistical distribution. 

A theory with the finite velocity of light, c, is called relativistic. The limit c — > oo 
leads to a nonrelativistic theory. A theory containing the Plank constant, h, is termed 
quantum. Taking the limit H — > 0, one passes to a classical theory. A theory dealing 
with probability distributions and, thus, containing the Boltzmann constant, k B , is 
named statistics. The limit k B — > reduces the latter to a particular kind of theory for 
which there is no a commonly accepted name. For such theories, we use below the word 
"mechanics", although this term is more general and is applicable not only to theories 
with k B = 0. So that the general name "mechanics" is to be considered not separately 
but in the combination with other terms characterizing a theory. The whole theoretical 
physics can be subdivided into 8 theories according as the principal parameters from 
the triplet {c, h, k B } are finite or some limits are taken. These theories are as follows: 

(1) Relativistic Quantum Statistics: {c, h, k B }; 

(2) Nonrelativistic Quantum Statistics: {c — > oo, ft, ks}; 

(3) Relativistic Classical Statistics: {c, h — ► 0, k B }; 

(4) Nonrelativistic Classical Statistics: {c — > oo, % — > 0, k B }; 

(5) Relativistic Quantum Mechanics: {c, h, k B — > 0}; 

(6) Nonrelativistic Quantum Mechanics: {c — > oo, H, k B — > 0}; 

(7) Relativistic Classical Mechanics: {c, H — > 0, k B — > 0}; 

(8) Nonrelativistic Classical Mechanics: {c — > oo, h — > 0, k B — > 0}. 

As is clear from this classification, Relativistic Quantum Statistics is the most gen- 
eral theory from which all other theories can be derived by means of appropriate limits. 
If not to subdivide the theories onto relativistic and nonrelativistic, then Quantum 
Statistics is the most general physical theory incapsulating all other theories. 
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